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SUMMARY 


Two different singularity methods have been utilized to 
calculate the potential flow past a three dimensional non-lifting 
body. Two separate FORTRAN computer programs have been developed to 
Implement these theoretical models, which will In the future allow 
Inclusion of the fuselage effect In a pair of existing subcrltlcal 
wing design computer programs. 

The first method uses higher order axial singularity 
distributions to model axi symmetric bodies of revolution in an 
either axial or Inclined uniform potential flow. Use of inset of 
the singularity line away from the body for blunt noses, and cosine- 
type element distributions have been applied to obtain the optimal 
results. Excellent agreement, to five significant figures, with the 
exact solution pressure coefficient value has been found for a 
series of ellipsoids at different angles of attack. Solutions 
obtained for other axi symmetric bodies compare well with available 
experi mental data . 

The second method utilizes distributions of singularities on 
the body surface, in the form of a discrete vortex lattice. This 
program is capable of modeling arbitrary three dimensional non- 
lifting bodies. Much effort has been devoted to finding the optimal 
method of calculating the tangential velocity on the body surface, 
extending techniques previously developed by other workers. Again, 



the best solution for ellipsoids at angles of attack ranging between 
0-30°, has been obtained using cosine spacing of the elements 
axially. 

For a simple axi symmetri c body of revolution without any slope 
discontinuity on the body surface, the first method offers a more 
accurate solution for less computational cost than the second 
method.* Therefore, while the surface singularity method has the 
advantage of greater geometry generality, the axial singularity 
method Is judged more suitable for modeling of fuselage effects in a 
preliminary aerodynamic design computer program. 


2 



CHAPTER 1 


INTRODUCTION 

Use of computers in aerodynamic theory has seen a trend in 
recent years toward the flow field analysis of complete aircraft 
configurations by a single program using a combination of boundary 
layer analysis and the full potential equation. Programs designed 
for this purpose do in fact exist, (refer here to references by 
Streett (1981) or by Mason (1977)), but most of them are extremely 
large and abound with subtleties often not evident to the user. In 
an ongoing research effort to include the body effect in a pair of 
existing wing design computer programs (Kuhlman and Shu, 1981; 
Kuhlman, 1983), it was decided to first develop separate fuselage 
programs as a preliminary step. Based on this spirit, the present 
work investigates the flow field around certain fuselage-shaped 
bodies using two theoretical formulations and their resulting 
computer programs. Either of these codes might possibly then be 
used in the above mentioned subcritical wing design programs to 
include the effects of a fuselage on the design. 

The problem under consideration in -this study is that of the 
uniform, incompressible, inviscid, potential flow past a non-lifting 
body. Bodies of interest may be either axi symmetric or arbitrary, 
although results to be presented herein have been limited to 
axi symmetric bodies. The governing differential equation in terms 
of the velocity potential $ is the well known Laplace's equation. 
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The surface of the body can, in general, be described by a function 

1* ^) ~ 0* The variables x^, x 2 , and x 3 are determined by 

the coordinate system used in the calculation. The mathematical 

statement of this problem is to find a solution to: 

V 2 * = 0 (1.1) 

subject to the boundary conditions 

grad $ • grad F = ^ • grad F = 0 (1.2) 

on the body, F(xj, x 2 , X3) = 0 
and 

grad $ = T at infinity (1.3) 

If the problem is rewritten in terms of the perturbation field 
<fr, due to the body, equations (1.1 to 1.3) become: 

V 2 * = 0 (1.4) 

(grad <j> + U^t) • grad F = 0 

on F(x-, ,x 2 ,x 3 ) = 0 (1.5) 

grad <j> -► 0 at infinity (1.6) 


Despite the fact that Laplace's equation is one of the simplest 



and best known of all partial differential equations, the number of 
known useful exact solutions Is quite small. The only exact 
solution of the uniform potential flow about a closed three- 
dimensional body is that for the general ellipsoid and its 
specializations. This is due to the difficulty in satisfying the 
nonhomogeneous boundary conditions and the requirement for the body 
surface to be a coordinate surface of one of the special orthogonal 
coordinate systems for which Laplace's equation can be separated 
into ordinary differential equations. As a result, use of numerical 
approximate solutions is usually considered in practical 
applications. Many different numerical schemes have been 
formulated. Some are exact in problem formulation and numerical 
methods then are used to obtain the approximate solution. Under 
such schemes, errors in the calculated solution can in principle be 
made as small as desired, by sufficiently refining the numerical 
calculations. In contrast, some schemes introduce analytical 
approximations into the problem formulation itself, and thus place 
limits on the accuracy that can be obtained regardless of the 
numerical procedures used. One typical example is the so called 
slender body theory where the body boundary condition is greatly 
simplified for high fineness ratio bodies. 

The present study is based on the first classification which is 
exact in formulation and numerical in solution. There are two main 
classes of numerical approaches; axial and surface singularity 
methods, which may be used in solving the potential flow past non- 
lifting bodies. From a practical point of view, these two methods 
each have their own characteristic properties and advantages. 
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Choice of method for a particular application is largely dependent 
upon the configuration of interest and how accurate the solution 
must be. In general, the axial singularity method has less 
flexibility than the surface singularity method, in that it can not 
handle certain axi symmetric bodies with discontinuous slopes, nor 
can it handle the entire class of non-raxi symmetric bodies. However, 
for simple and axi symmetric body shapes, the axial singularity 
method is easier and more economical to use than the surface 
singularity method. In contrast, the surface singularity method is 
in principle more accurate and has fewer geometrical restrictions 
when compared to the axial singularity method. It is, then, much 
more sophisticated and powerful, but requires a substantially 
greater number of numerical calculations and more computer 
storage. Other solution methods applicable to solving Laplace's 
equation include finite difference methods (David and Geppson, 

1973), but these methods have not been utilized in the present 
study. This is because the singularity methods based upon Green's 
theorem significantly reduce the cost of obtaining a numerical 
solution for <j> to a certain accuracy, largely as a result of 
reduction of the problem to finding a suitable singularity 
distribution on the body surface, rather than requiring that a 
solution for <j) be determined in the entire region exterior to the 
body. The development of the present singularity methods have been 
summarized in Chapters 2 and 3 of this thesis, for the axial and 
surface singularity methods; respectively. Also given in these 
chapters are brief literature surveys. 

In light of these considerations, a FORTRAN computer program 
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based on the axial singularity method has been implemented in the 
present work for the simple, ax i symmetric fuselage. Much effort has 
been devoted to obtain the most accurate numerical solution for both 
axisymmetric flow and inclined flow. Higher order axial singularity 
distributions have been formulated for both axial flow and cross 
flow. Pressure distribution results for ellipsoidal bodies as 
presented in Chapter 2, show excellent agreement with the exact 
solution everywhere on the body at arbitrary angle of attack 
providing that separation does not exist on the body. In Chapter 4, 
additional’ results from the axial singularity method, for more 
complicated geometries, have been compared with results from a 
surface singularity program described in Chapter 3 and with 
experimental data. 

Although the surface panel method is not considered an 
efficient model to evaluate fuselage effects in the preliminary wing 
design program, because large computer memory and long computational 
time are needed to perform the calculation, it does offer a reliable 
solution which can serve as a bench mark for those solutions 
obtained from the axial singularity method. Chapter 3 gives a 
detailed outline of the formulation and accuracy of the surface 
panel technique utilized in the present study. A second FORTRAN 
computer program has been implemented as part of the present work to 
perform numerical calculations, for ellipsoidal bodies in 
axisymmetric flow or at angle of attack. These results have been 
presented in Chapter 3, while results for more complex geometries 
have been presented in Chapter 4. 
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LIST OF SYMBOLS 


a 

A 

b 

Cp 

d 

e 

i 

L 

Ka 

Kc 

M 


n 

N 

p ij’ q 1j’ r ij 

q 


Inset distance or semi-major axis of ellipsoid 
area of vortex lattice panel 
semi -mi nor axis of ellipsoid 
pressure coefficient 

distance from control point to vortex filament 
eccentricity of ellipsoid 

body length or vortex filament length or panel chord 
length 

total body 1 ength 

axial coefficient of virtual mass in equn. (2.66) 

cross coefficient of virtual mass in equn. (2.66) 

order of singularity strength polynomial function or 
Mach number . 

unit normal vector 

total number of axial singularity panels or total 
number of vortex lattice panels 

vortex lattice influence coefficients, see equn. (3.3) 
source strength 


^ radial coordinate or magnitude of position vector 

R nose radius of curvature or shoulder radius of body 

s bounding surface of a region 

u axial velocity component 

u x’ u r ,u 8 axial, radial and circumferential velocity components 

U. free-stream velocity 
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V 


w 

W 


x,y,z 

x,r,8 

a 

3 

r 

0 

0 

U 

5 

IT 

P 

T v 

$ 

* 

SUBSCRIPTS 

c 

1 

j 

k 


total tangential velocity or induced velocity due to a 
vortex filament 

total tangential velocity 

radial velocity component 

self-induced velocity due to a closed loop with unit 
strength of circulation 

cartesian coordinates of points 

cylindrical coordinates of points 

angle of attack 

inclination angle of control point with x-axis 

strength of circulation 

circumferential velocity component 

circumferential coordinate or meridian angle or 
inclination angle of panel with x-axis 

doublet strength or curve! inear coordinate variable 

axial integration variable 

ratio of circumference to diameter of circle 

radial distance of control point to body axis 

tangential velocity due to local vortex lattice panel 

velocity potenial 

disturbance velocity potential 


index of panel control point or centroid of area 
index of control point for axial singularity method 
index of influencing panel 
summation index 
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index of minimum induced velocity point 

index of axial and radial components for axial flow 

index of axial, radial and circumferential components 
for crossflow K 


ABBREVIATIONS 


CS 

constant source 

CD 

constant doublet 

LVS 

linearly varying source 

LVD 

linearly varying doublet 

CPE 

exact value of pressure coefficient 

CPN 

numerical value of pressure coefficient 

SR 

slenderness ratio 
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CHAPTER 2 


AXIALLY DISTRIBUTED SINGULARITY METHOD FOR 
AXISYMMETRIC BODY OF REVOLUTION 

2.1 Literature Survey 

The idea of utilizing a distribution of singularities interior 

to the body surface to represent the potential flow past a body of 

revolution was first introduced by Von Karman (1930), who considered 

axi symmetric shapes in axi symmetric flow and represented them by a 

piecewise constant source distribution along the axis of symmetry. 

Since the original work of Von Karman, very little information on 

this method has been generated, until recently when Oberkampf and 

Watson (1974) discovered that this method produces a system of 

linear equations which is in general ill-conditioned, and requires 

very high computational accuracy (25 significant figures) in the 

construction of the influence coefficient matrix. They concluded 

that the method does not always produce reliable solutions for the 

flow around a specified body and is very sensitive not only to the 

shape of the contour but also to the number of elements used to 

generate the body. Karamcheti (1966) states that the body should be 

slender and should not have any discontinuities in the slope of the 

meridian line. Zedan and Dalton (1978) extended the method by 

employing piecewise linearly varying source distributions to improve 

the computational accuracy. The solutions obtained therein are more 

accurate as compared to the constant strength element methods. One 
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of the most important features of the linearly varying source method 
is its ability to deal with bodies having an inflection point in the 
meridian contour. They again concluded that their new method cannot 
handle body shapes with sudden changes in the slope of the meridian 
line. More recently, Zedan and Dalton (1980) have used a polynomial 
function of arbitrary degree to represent the variation of the 
intensity of the source distribution over each element. The effects 
of the order of the distribution, the number of elements, the 
normalization of the body coordinates, the fineness ratio and the 
geometry of the profile on the performance of the method have been 
studied in detail. They concluded that with appropriate choice of 
these parameters, this approach can be as accurate as the surface 
singularity approach even for simple axi synanetri c bodies with 
inflection points. 

In view of the simplicity of the axial source distribution 
methods as compared to the surface singularity distribution methods, 
and the tremendous saving in computer memory storage and the amount 
of numerical calculation necessary, it was considered justifiable to 
apply this method as an initial step in considering the effect of an 
axi symmetric fuselage upon the optimal camber surface obtained from 
an existing aerodynamic wing-design computer program, provided the 
accuracy can be improved and the limitations can be reduced for the 
former class of methods. 

The present study is an extension of Zedan and Dalton's work. 

An axial singularity distribution is placed along the centerline of 
the body of revolution, where this singularity distribution is 
assumed to have either a piecewise constant or a piecewise linearly 
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varying strength over a series of discrete segments. Expressions 
for velocity potential and induced velocities have been developed 
for source (identical to those of Zedan and Dalton, 1980) and 
doublet (as suggested in Karamcheti, 1966) singularity 
distributions. In addition to this, the expressions have been 
extended for piecewise constant, or linearly varying crossflow 
doublet singularities, which allows simulation of an axisymmetric 
fuselage at angle of attack. 

A FORTRAN computer program has been developed to Implement this 
theoretical model. This program is currently available for use in 
the optimal wing-design computer programs (Kuhlman and Shu, 1981; 
Kuhlman, 1983). Future plans are to add the current axisymmetric 
body model in the wing-design code. In this chapter the development 
of the theoretical model will be summarized, and the effects of the 
type of singularity, the order of the distribution, the number of 
elements, the geometry of the profile, as well as some paneling 
techniques, on the performance of the method will be shown by using 
a simple ellipsoidal test case. Results have been obtained for a 
series of ellipsoidal bodies of revolution of varying slenderness 
ratio and angle of attack, and these results have been compared with 
an exact solution summarized in Wang (1970). Results for more 
complicated axisymmetric bodies will be presented in Chapter 4. 

2.2 Formulation and Equations of the Problem 
2.2.1 Formulation of the Problem 

Consider the steady, acyclic, potential flow past an 
axisymmetric body of revolution that is at an angle of attack to the 
direction of the undisturbed stream (see Sketch 2.1). 
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Sketch 2.1 Illustrating flow with angle of 

attack past the body of revolution. 

In analyzing the flow past a body of revolution, it is 
convenient to use cylindrical coordinates (r,8,x) as shown in Sketch 
2.1. The surface of the body is described by an equation of the 
form r * R(x). This may also be expressed as F(r,x) = r - R(x) = 0. 
Equations (1.4 to 1.6) then take the following form: 



( 2 . 1 ) 


(U^sinasine + U p ) + (U^coso + U x ) |£ = 0 
on F(r,x) = 0 


V* 0 at infinity 


where U = |4. , u =* 1 i_£ n - 

r 3r e r 36 » u x lx 


Due to the linearity of the problem (governing 


( 2 . 2 ) 

(2.3) 


equation and boundary 



conditions), the solution of Equations (2.1 to 2.3) may be expressed 
as the superposition of the solutions of two separate problems: (1) 

the problem representing the axial flow at speed U^cosa past the 
body of revolution, and (2) the problem representing the cross flow 
at speed U^sino past the same body of revolution (See Sketches 2.2 
and 2.3). 


Z 



Sketch 2.2 Illustrating axial flow past the body 
of revolution 


Z 



I I I U^sina 

Sketch 2.3 Illustrating cross flow past the body 
of revolution 
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These two mathematical problems may be stated in the form of: 
Axisymmetric Flow Past a Bocty of Revolution: ^ = ^(r.x) 

= 0 (2.4) 

U rl (R,x) - U xl (R,x) = U^coso -gB. (2.5) 

* 0 at infinity (2.6) 

Cross Flow Past a Body of Revolution: <p 2 = ^(r.e.x) 


? 2 *2 s 0 (2.7) 

" U X £(R»x) = -U^sinasine (2.8) 

7(J>2 0 at infinity (2.9) 

Then one is able to use line singularity distribution 
techniques to solve these two problems and obtain the total solution 
by superposition of 4>-j and <p 2 . For the axisymmetric flow part, one 
might use sources or doublets along the body axis. The singularity 
distribution is divided into finite elements. The strength 
variation over each element can be described by a polynomial 
function of arbitrary degree. In this study, we consider four types 
of representations, which are: constant source, linearly varying 

source, constant doublet, and linearly varying doublet. For 
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convenience to extend the method to higher order singularity 
expressions, the general formulas have been listed for reference. 

The cross flow part, as suggested by Karamcheti (1966), has been 
represented by a piecewise constant doublet or linearly varying 
doublet distribution along the axis of the body, with the axes of 
the doublets being normal to the body axis; the doublet axes oppose 
the direction of the undisturbed crossflow stream at infinity (see 
Sketch 2.4). After solving these two mathematical problems, one can 
obtain the pressure coefficient by using Bernoulli's equation 

V 2 

C 0 * 1 - ^ (2.10) 

P U c 
00 



Sketch 2.4 Cross flow past a body of revolution 

by combination of a doublet distribution 
and a uniform stream 


Although one can superpose the <j> and one cannot do the same for 
the pressure coefficient since it is a nonlinear function of 
Here V is local total velocity, in the following form: 
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V 2 = (U^cosa + U xl + U x2 ) 2 + (U^sinasinQ + U rl + U r2 ) 2 


+ (U sinacose + UJ‘ 

OO 0 


( 2 . 11 ) 


and 


3 ( j> ] 

*xl ~ Jx 


3<j>, 

I = _ i 
x2 3X 


* U rl = W 


d<pr 

I = — 
rZ 3r 


U 


_ 2^2 

e r 30 


(2.12) 


2.2.2 Equations for an Element 

The following three techniques have been used to evaluate the 
integrals which will appear in the induced velocity expressions: 

(1) Integration by analytic geometry relationships, (2) integral 
tables, and (3) the MACSYMA symbolic manipulation language (Bogen et 
al., 1975). All three methods have been utilized in the present 
work, with the majority of the integrations having been performed by 
MACSYMA. 

In the following six different element equations, we define the 
symbolic expressions: 



c K-1 


(2.13) 


*2 ^K-l 

>K * N d5 


(2.14) 


where U^tg) are the X coordinates of the endpoints of the 
element. Then the MACSYMA symbolic manipulation language has been 
utilized to obtain the integral expressions as listed in Appendix A. 
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2.2.2. 1 Constant source distribution. Consider first 

the j** 1 element of an axial singularity distribution having a 

constant Intensity q. as shown in Sketch 2.5. The velocity 

J 

potential function and the velocity components at 
the 1 th point (x. ,r f ) are given by 



w(Xj ,r.. ) 


L 1 C(x i -?r + r. d T 


(2.17) 


Here q. is constant on element j. 

J 

Using tables of integrals (Gradshteyn and Ryzhik, 1980), u and 
w are evaluated in closed form as 






1 


w(x i’ r i> b -37F7< 


. q jl , £ 2 ’ x i 


h ' x i 


1 V* 2 Z " 2 V2 +x i + r i ' 2x i £ l + x i 2 + r i 2 


(2.19) 


In MACSYMA language expressions (Appendix A), the constant source 
velocity may be expressed as 


u(Xj ,r.j ) 





( 2 . 20 ) 


w(x 1* r 1 J 


r 1 H 1 


(2.21) 


2. 2. 2. 2 Linearly varying source distribution. If the j th element 
(see Sketch 2.5) Is assumed to have a linearly varying source 
distribution, q^ = a^ + a^S* the velocity potential function and 
the velocity components at the 1 th point (x f ,r^) are given by 

r 


U cosa 
00 


source strength q. 

J 

Sketch 2.5 Illustrating axial singularity panelling scheme 



20 


.r, ) • - k l [ z aj2 s d< 

Vlx ,.-!) 2 7 r , 2 


1 *”? ^ll +a i? 5)(x.-?) 


1 C (x 1 -?) 2 + r. 2 ] 3/2 


r. z ? a., + a. ? c 

w(x i’ r 1 ) = ■?? V — ~ 2 J2 2 3/2 d5 
11 4,r *1 C(x r c) 2 + r^ 2 ] 3 ' 2 


( 2 . 22 ) 


(2.23) 


(2.24) 


Equations (2.23 and 2.24) can be rewritten as: 

(2.25) 

(2.26) 


u( V r 1> = -Jf C V X ! H r H 2 ) + a j2 (x 1 H 2- H 3 )] 
w( V r 1> = W C a jl H 1 * a j2 H 2^ 


At this point, to extend the scheme to an arbitrary order 
polynomial function, simply let qj =[?£] a^s*" 1 , substitute back 
into equations (2.16 and 2.17) and evaluate the integrals to 
obtain the following general expressions: 


u(x 1- r 1> - -W !S1 a JK (x 1 H K - W (2-27) 

jr* # 

w(x i* r i) = W k£] a jK H K ( 2 - 28 ) 

where a^ are constants. 
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2. 2. 2. 3 Constant doublet distribution. Next, the source 
singularity is replaced by a doublet singularity with axes directed 
into the undisturbed stream (see Sketch 2.6). That is, the doublet 
strength Is written as = -Ujf. The velocity potential function 
and the velocity components at the i th point (x^ ,r.j) due to 
the element are given by: 


♦^ x 1 » r j) ~ I o 5 TT7T 

' 1 ” *1 C(x r ?) 2 + rff 72 


dS 


(2.29) 




3 r l 2 "j (x r 5) ' 


37 ? 


dC 


u /v . ) = _ 3r i / 2 y j (x i" 5) 

1 4T h [(x,-S) 2 + r. 


t^tz d5 


(2.30) 


(2.31) 


Equations (2.30 and 2.31) may be evaluated for a piecewise 
constant doublet distribution in the following forms: 


u - 3u . 

= — J u _ J ty, 


U(x i* r i ) = # H 1 - (x i P 1 " 2 V 2 + P 3 ) (2.32) 


3r,u, 


w(x 1- r 1> = — 3GT <x ( p l - p 2> 


(2.33) 
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U COSa 
00 


r 



Sketch 2.6 Illustrating axial doublet singularity 
panelling scheme 


ith 


2. 2. 2. 4 Linearly varying doublet distribution. 

If the doublet strength Is linearly varying over the j 1 
element, -Mjt * -(a^ + a j 2 5 ^* the potential function and 

the velocity components at the 1 th point (x^.r^) due to this panel 
are given by: 


♦(xj.rj) 


u(x i ,r < ) 


i ‘s 

' l-t 777 r\2.\ - 2 t < 3/2 


1 C(x r 5) 6 + rfr 


d? 


(2.34) 


1 H 


a jl + a j 2 5 


tlx ,-?) 2 7 d5 ' * ; ‘1 


3 f l 2 (a jl* a j2 g)(x 1~ g) ' 


C (x 1 -c) + rp 


7372" 


d5 


(2.35) 


w(x.j,r.|) 


3r 1 r l 2 (a J1* a J2 g)(x r e) 


tr 


f 


*1 C(x 1 -6) 2 + r^] 


27377 


d5 


(2.36) 


By evaluating the above expressions, again using MACSYMA, one 23 



obtains: 


u(x i’ r i ) = ir Ca jl (H l " 3x f 2p i + 6x^2 - 3P 3 ) 

+ a j2 (H 2 - 3x 1 2p 2 + 6x i P 3 - 3P 4^ (2.37) 

3 r. 

w(x 1’ r i ) = JT Ca jl (P 2 ' x iV + a j2 (P 3 “ x j P 2 (2.38) 

^ U ^+1 order doublet distribution is necessary, use 

~ u jt = “ K ^ a jK 5 ^ t0 re P lace Uj in equations (2.30 and 2.31), 
and then the general expressions for ($,u,w) are: 


, M+l 

* (x i> r i ) * • a jK (H K+l - W 


(2.39) 


u(x i* r i> = *r K f, V H K - 3 x i Zp K + 6 Vw 3 W (2 - 40 > 


M+l 


3r. M+l 

w(x ) ’ r i 1 = TT K 2 , a jK^K+l ” x ) 2 K> (2-41 ) 

2. 2. 2. 5 Constant cross flow doublet distribution. Use of a 
piecewise constant doublet distribution along the body axis to 
simulate cross flow has been illustrated in Sketch 2.7. 

The j element (from x * z.j to ^) is a constant strength 
doublet distribution with doublet axes pointed toward the on-coming 
flow. The velocity potential function * at a field point "i" due to 
this element can be expressed as: 
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Sketch 2.7 Illustrating axial doublet singularity 
‘panelling scheme for cross flow 


^ij( x i ,r i , 9 i) 


sine.. 

TT" 


U 


2 u j r 1 


d5 


(2.42) 


The reason that is involved in the velocity potential 
function expression is due to the three dimensional behavior of the 
cross flow. This 0^ dependence will be kept in order to make the 
expressions derived complete. However, later on, it will be shown 
that the doublet strength distribution gotten from solving the 
linear equations resulting from the boundary conditions is 
independent of the e location of control points. In this sense, we 
can eliminate tne e^ dependence simply by locating control points on 
the top meridian line (e^ = 90°). However, if it is desired to cal- 
culate the pressure coefficients at points not located on the top or 
bottom meridian lines, the 0 effect must be included. 

By differentiating three velocity expressions result: 

* J 
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u(x 


3sine. ju u.r.(x.-c) • 

ds (2 ' 431 


sine. i ? 

w{x i’ r i’ 9 i> 


"j 


1 [(»,-«)* + T^F 2 


d5 


3sine. £ 0 

\ t L 


“j r i 


**“ ! h [(x r5 ) 2 1 7^i d5 


(2.44) 


cose. 

0 ^ x 1» r 1* e f) = 1 r 2 


4tt J Jt. 


W J 


1 C(x r e) Z + r. 2 ] 3/2 


d5 (2.45) 


If the procedure described in section 1.2.2 is followed, the 
expressions for velocities will be obtained as follows: 


u.sine. . 

u ( x 1 ,r 1 ,e f ) .-J’ (sl „ 3 6 j>j+ i. sj „ 3 6nj) 


4irr i - 


(2.46) 


w(x 


u - sine 


i ,r f’V = i^rr C 2 (cos 6 i.j*i - C0SB U> - (cos 3 s i,j + i - cos 3 ei>j ): 


(2.47) 


i^cose . 

(x^r^e.) _ +i 


(2.48) 


where the 0 1 -j +1 » 3^ j have been defined in Sketch 2.7. 

In the computer program, a subroutine based on these 
expressions has been implemented. If the MACSYMA language is used 
to perform the above integrations, the following forms result for 
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the constant cross flow doublet: 


3u .r .sine . 

“ ( *1" r 1- 8 1> T7 (x l p l - p 2> 


. 2 
u.sine. 3u.r. sine. 

w(x i- r i- 9 i> • ■ " V" H i hi p i 


(2.49) 

(2.50) 


0 (x^rj.Oj) 


u .cose . 

J ' u 

— *r- H i 


(2.51) 


2. 2. 2. 6 Linearly varying cross flow doublet distribution. If the 
doublet strength is linearly varying over the j element, 

-u,J = -<a„ ♦ a j 2 ?)j , the velocity potential function and the 
velocity components at the i th point (x. .r^ ,0^ ) due to the j th panel 
are given by: 


sine. 


r*2 (a jl + a j2«> r 1 

\ tCTTTW ds 


(2.52) 


u(x l .,r.,e i .) = - 


3sin8 i *2 (a jl +a j2« )( V 5)r t 


"4F 


A 1 C(x t -€) 7 + r. 2 ] 5/Z 


d£ (2.53) 


wfx^.r^e^ 


sine.. 

4tt 


*2 a ,11 * a i2 ; 


d5 


3s(n6 i f *2 <a jl'*' a j2 5,r l 

’ ! h [(x r5 )5 ♦ r « 


cose, i. 


0 (x^r^e.) = — 2 


a jl + a j2 5 . 


1 C(x.-c)‘ + r. c T 


(2.54) 

(2.55) 
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By using the MACSYMA language, one obtains: 


3r .sine. 

u(x C r fV * H Ca jl < x f p 1 - p 2 > * a j2^ x 1 P 2 * P 3 )] (2 ' S6) 

si ne . 9 

w(x 1> r f e 1 ) ■‘-Sr-C a jl< H l - 3r, z P,) + aj2 (H 2 - 3r,. 2 P 2 )] (2.57) 
cose . 

® <*, .-><>,) = ^-L [ ajl H, + a J2 H 2 ] (2.58) 

The same method applied to arbitrary higher order doublet 

M+l j, . 

distributions, using -u.j = - z a..5 j to replace y. in 

J K=1 J 

equations (2.43 to 2.45) results in: 

3r.sine. M+l 

“‘WV 3ir a jK (x 1 P K - P K+1 ) < 2 -5 p ) 

sine. M+l 

w(x f .ri,0 f ) a JK (H K - 3r f 2 P K ) (2.60) 

cose. M+l 

© (x.,r.,e.) = s^— ^ a jK H K (2.61) 

2.2.3 Equations for Determining the Strength of Singularities 
The expressions for induced velocity components at 
the i th control point due to the j th element distributed along the 
body axis (from x = to l^) have been derived in section 2.2.2. 
When the contributions from individual elements are added together, 
the total induced velocity at the i control point can be expressed 
in terms of the unknown strengths of the singularities at each 
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element. These unknown strengths can be specified as arbitrary 
order polynomial functions. In the present studies, constant (order 
M = 0) and linearly varying (order M = 1) cases have been 
considered. When these expressions are substituted back into 
equation (2.5) or equation (2.8) to satisfy the no penetration 
requirement on the body surface, a system of linear algebraic 
equations for the unknown polynomial coefficients is obtained, which 
is then solved for the unknown strengths. Once the strength 
distribution is known, the velocity and pressure coefficient can be 
calculated at any field point outside the body, as well as on the 
body surface. 

Since the solution procedures are essentially the same in 
either the axial flow or cross flow calculations, therefore, only 
the source singularity distribution in axial uniform flow is 
presented as an example. There is a slight difference in using the 
linearly varying doublet as the singularity in either the axial or 
cross flow calculation. This difference will be described 
subsequently. 

The induced velocity components at the i th control point due to 
N elements distributed along the body axis (see Sketch 2.8) are 
given as follows: 


1 N 

u( V r i> ' 

J * 


M+l 

a jK (x i H K " ^+1^ 


(2.62) 


r . N 

w(x i’ r i } = W * 

J ' 


M+l 

C K=l 


(2.63) 
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where M is the degree of the polynomial function, and a^ is the 
unknown coefficient of the singularity distribution element "j". 
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Sketch 2.3 Illustrating axial source singularity 

panel distribution for equations (2.62 & 2.63) 

If the constant source is chosen as the singularity type(M = 

0), there are N unknown coefficients and a N*N linear system results 
from applying the no penetration boundary condition at the selected 
H control points on the bo<Jy surface. Gauss elimination methods 
have been employed herein to solve for the unknowns. If M = 1, 
linearly varying source. Is chosen as the singularity, there are two 
unknowns for each element. However, continuity of the singularity 
distribution at the junction point between two neighboring elements 
provides another N-l equations. This continuity requirement Is 
valid so long as the bocfy slope and curvature are continuous. The 
general form of these N-l continuity equations Is given as: 

a j,l + a j,2 A j * a j + l,1 + a j+l ,2*j J * 1.2 . • . N-l (2.64) 

where Is the x coordinate of the junction point of elements "j" 
and "j+l\ 

The last equation required is a closure condition, due to the 
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fact that the net efflux of the source distribution should be zero 
for a closed body, and is given by 


(a j,l L j + 7 a j,2 L j 2) = ° 


(2.65) 


.th 


where Lj is the length of the j element. Hence, the 2N unknowns 
may be determined from the 2N equations developed above. If higher 
order singularities are used (M > 2), then corresponding higher 
order derivatives of the source strength should be continuous, at 
least for bodies which are sufficiently smooth. For example, for 
the second order source (M 3 2), continuity of the slope of the 
source strength should be enforced at the juncture between two 
neighboring elements. 

As mentioned above, when the linearly varying doublet 
singularity is used to represent the body, the last closure 
condition is no longer applicable. The form of this last condition 
for a general body of revolution is not known. However, in the 
present studies, a condition of zero strength at the starting and 
ending points of the singularity distribution has been used to 
perform the calculations. Imposition of some more fundamental 
condition on the linear doublet distribution may improve the 
solution accuracy for this singularity type. 


2.3 Analytical Solution for Ellipsoidal Body 
The potential flow past an ellipsoidal body was solved exactly 
by Lamb (1932). Wang (1970) gave an explicit expression for the 
surface velocity potential function * in terms of two surface 
coordinates u and 8 (see Sketch 2.9). Here, (-1 < u < 1) is 
constant along the parallels and 9 is constant along the 


31 



meridians. It should be understood that no restrictions on angle of 
attack or slenderness ratio was imposed in the solution. By 
differentiating the velocity potential, *, and utilizing the 
Bernoulli's equation, a closed form of the surface pressure 
coefficient, C p , is obtained as: 



Sketch 2.9 Coordinates scheme for equation (2.66) 


, H * K a )(cosa)(l - u 2 ) 1/2 + (|)(l + K c )(s1na)iiCose 


c p * 1 -<c 


+ CO + K c )s1na sine] 2 ) 


(1 - eV) 


( 2 . 66 ) 


2 

e 3 d , the eccentricity 

di 

^ log - 1 

K a * 7*"1 I ^ T+eT ' the ax * aT coefficient of virtual mass 
C WT * 2? ’<>3 TT? 

K c a T+2K » the cr °s* coefficient of virtual mass 
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A similar, explicit expression for Cp has been given by Cebeci , 
Kattab and Stewartson (1980). Another expression for Cp for an 
ellipsoid at angle of attack has been given in Schlichting and 
Truckenbrodt (1979), expressed in cylindrical coordinates (see 
Sketch 2.10). The pressure coefficient, Cp, was expressed in the 
following form: 



Sketch 2.10 Coordinates scheme for equation (2.67) 


c p- [1 


where 


1 - (*) 2 * 
A 2 


l-Cl - \](|) 
a^ a 


/ x 2 
VI -(-J) cose 
T \ 2 o (2.67) 
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S 2 
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.2 1/2 
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) ] - 


b 2 1/2 

(i - K) > 

a 


It is easy to show that both expressions are identical for flow 
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without angle of attack (a = 0°). However, an appreciable error 
results from equation (2.67) for the case of the flow at an angle of 
attack (see Table 1 and Figure 1). Numerical experiments show; 
at * a 5°, the discrepancy from the exact solution (Lamb) is about 
7%; while at o = 10°, more than 30$ discrepancy results, which is 
totally unacceptable in any kind of flow calculation. Schlichting 
and Truckenbrodt indicated that this equation is only good for flow 
at small angles of attack, but such a large discrepancy is 
unexpected. It was found that the small a equation (2.67) can be 
directly derived from the exact equation (2.66) by assuming a is 
small, so that cosa ~ 1 and sina ~ a, and neglecting all terms of 
order a . Therefore, while equation (2.67) is the correct 
expression for Cp for small a, use of equation (2.66) is recommended 
instead of (2.67) for all values of a. In the comparisons between 
the present numerical model and the analytical solution to be 
presented in the following section, only the full expression given 
by equation (2.66) has been utilized. 

2.4 Accuracy of Method 

It has been found that the solution accuracy of the axial 
singularity method is largely dependent upon the following four 
factors; (1) singularity type, (2) number of discrete elements, (3) 
inset panelling technique, and (4) panelling scheme. In this work, 
a systematic evaluation has been conducted to determine the error 
sensitivities of the solution for an ellipsoidal body to these 
factors. The motivation to do these studies is basically due to the 
fact that decisions on these variables must be made a priori for any 
particular problem, and will usually be influenced by the problem 




and the information desired from the solution. It is expected that 
the present study might provide some guidelines for more general 
bodies. 

The test geometry has been chosen to be an ellipsoid of SR=5, 
because the exact solution previously described already is known, 
and it is a suitable geometry to show the error sensitivities for 
all four factors mentioned above. Similar studies of SR=2 and SR=10 
ellipsoids have yielded solution of accuracies comparable to those 
presented herein for SR=5. All comparisons have been made using 
root mean square error for comparison. This error calculation has 
been based on pressure coefficients calculated at 40 points equally 
distributed along a meridian line on the body surface. The root 
mean square error is defined as 



(CPE. - CPN f ) 2 
TO 


( 2 . 68 ) 


where CPE. is the exact pressure coefficient at the i tfl point 

1L 

and CPN^ is the numerical pressure coefficient at the i point. 

The same fractional locations of the 40 points have been used, 
independent of the control point locations. 


2.4.1 Test 1: Axi symmetric Flow 

A set of results for 19 typical cases for the SR=5 ellipsoid 


at a = 0° are listed in Table 2. Through comparisons of the error 
magnitudes, one can easily find the individual contributions to 
solution accuracy due to one factor by holding the others fixed. 
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2. 4. 1.1 Element number. By comparing configurations (1 and 3), it 
is seen that the solution is improved by utilizing more elements to 
represent the body. Also, the pair of configurations (2 and 9) 
reinforce more dramatically the same conclusion. In fact, 
increasing the number of elements generally improves the solution 
accuracy for any numerical calculation procedure. In the present 
studies, 20 elements are considered sufficient to represent the 
ellipsoid. 

2. 4. 1.2 Singularity type. By comparing configurations (5, 10, 12 
and 13), it may be seen that the constant doublet is the least 
accurate among the four types of singularity, being inadequate to 
represent the source like nature of flow near the stagnation 
point. Therefore, it has not been considered in further 
application. Next, the constant source distribution has the same 
accuracy as the linearly varying doublet. Finally, the linearly 
varying source yields the best solution, in the test case, for a 
fixed number of elements, and is recommended for application to 
other more complicated geometries. 

The singularity strength and Cp distribution along the body 
axis for these four configurations have been shown in Figures 2 and 
3. Note that the source strength predicted by slender body theory 
has been shown in Figure 2 for comparison. Since the numerical 
solution error in Figure 3 is almost indistinguishable. Table 3 has 
been prepared to display the deviations of these four solutions from 
the exact solution. Note that the error in Cp for the linearly 

_5 

varying source distribution is on the order of 10 over much of the 
ellipsoid, with larger errors near the stagnation points. 
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As a fairer comparison of the relative accuracies of the 
constant and linear source singularities for approximately equal 
computational cost, configuration pairs (13, 17) and (14, 18) may be 
compared. These test pairs have equal matrix sizes and thus should 
have nearly equal run times and costs. It is seen both for equal 
sized element (configurations 13, 17) and cosine spacing 
(configurations 14, 18) use of 40 constant source singularities 
leads to better accuracy than 20 linear sources. However, 
comparison of pairs (5, 19) displays an opposite trend for a sparser 
element distribution. That is, ten equal sized linear sources yield 
a more accurate solution than 20 constant source singularities. 

Based on these limited studies it is apparent that use of higher 
order singularities may or may not result in reduced computational 
cost. 


2. 4. 1.3 Panelling scheme. The semi-circle cosine type panel 
distribution has been used previously by many authors in choosing 
the chordwise and spanwise locations of vortex points as well as 
control points (see, for example, Lan, 1974). This panelling scheme 
is based on the assumption that better results might be obtained if 
finer elements are used in regions of rapid variation of sectional 
properties. Of course, qualitative properties of the section must 
be known a priori. It has been proven effective and accurate 
through numerical experiments in thin wing theory. Therefore, it is 
believed to be reasonable to apply this technique in the present 
body calculation. Two comparisons have been made; first, between 
configurations 5 and 9, and second, configurations 13 and 14 of 
Table 2. Both comparisons give the positive answer, that a 
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tremendous improvement in solution accuracy is obtained from using 
cosine spacing instead of an equal spacing scheme. Hence, use of 
cosine spacing is strongly recommended, wherever applicable, in 
other flow calculations. 

2. 4. 1.4 Inset panelling technique. Zedan and Dalton (1980) 
concluded that use of the submerged panelling technique would lead 
to a more accurate solution for bodies with either a blunt nose or 
tail. In the present studies, it has been found that the optimal 
inset distance for the ellipsoid of SR=5 is about 1% of the body 
total length (configurations 3 to 7). The root mean square error 
versus submerged distance is shown in Figure 4. 

For ellipsoids of different slenderness ratios, the optimal 

inset distance must be varied with the magnitude of the local radius 

of the nose and tail. This ellipsoid nose radius is found to be 

b/SR, where b is the length of the minor axis of the ellipsoid. 

Numerical experiments have shown that 6 % of the body length is the 

optimal inset value for a SR=2 ellipsoid, and 0.25% of the body 

length is optimal for a SR=10 ellipsoid. Table 4 (a), (b) display 

the relation between the inset value and the associated root mean 

square error. All test configurations have been run using 20 cosine 

space linearly varying source elements. 

« 

It has been discovered that a nearly constant relationship 
exists for the optimal inset value for other slenderness ratio 
ellipsoids, as shown the following brief table, prepared for 20 
cosine spaced elements. Here R is the nose radius of curvature. 
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TABLE 2.1 Optimal inset value for various SR ellipsoid 


optimal inset of singularity line 


SR a_ 

10 10 

5 5 1 0.2 1 

2 2 1 0.5 6 


% of a 

% of b 

(% of b)/R 

0.5 

5. 

cn 

O 

• 

2. 

10. 

50. 

12. 

24. 

48. 


b_ R_ % of L 
1 0.1 0.25 


Notice that the optimal inset value is a constant fraction of 
the ellipsoid nose radius. Therefore, the optimal inset value for 
all ellipsoidal bodies should be governed by the following 
expression: 

optimal inset value in terms of _ 25 (2.69) 

the body total length (%) ” "^7 

Similar results were found using 20 equal sized elements, but the 
constant was changed slightly so an effect of the length of the 
element nearest the nose could be seen. 

Table 4 (c) shows the root mean square Cp error for 
axi symmetric flow past a series of ellipsoids using the suggested 
inset value calculated from equation (2.69). Undoubtedly, this 
equation offers a very accurate Cp solution, and it's use is 
recommended in further calculations for blunt nosed ellipsoidal 
bodies. Figure 5 displays the Cp distributions for axi symmetric 
flow past SR=2, 5 and 10 ellipsoids using the optimal inset 
panelling scheme. 
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2. 4. 1.5 Optimal solution. By combining all of the above favorable 
factors together, an optimal representation for axi symmetric flow 
past a SR=5 ellipsoid may be determined as 20 discrete elements of 
cosine spaced linearly varying source distribution with a 1% inset 
panelling scheme. This configuration does give the minimum error in 
Cp for 20 elements of 0.000063. Figures 6 and 7 display the 
pressure coefficient and the error distributions. Note that the 
numerical solution matches the exact solution very well on the whole 
body surface, with the largest errors occurring near the nose and 
tail. Also shown, for completeness in Figure 6, is the Cp 
distribution calculated by slender body theory using the source 
distribution shown in Figure 2. Slender body theory is in error 
nominally by 15% over the middle 80% of the body. Similar 
calculations at SR=10 show approximately 4% error. 

2.4.2 Test 2: Inclined Flow at 5° Angle of Attack 

A procedure similar to that performed for a = 0° above, is now 
repeated for a = 5° to investigate the error sensitivities of 
panelling and singularity variables in the cross flow calculation. 
Since the inclined flow can be found as the superposition of a pure 
axial flow and the corresponding cross flow, one can fix the 
parameters of the axial flow part of the calculation, and allow only 
the cross flow calculation singularity and panelling to vary. A 
choice of 20 discrete elements, for a cosine spaced constant source 
distribution with 1% inset panelling will allow sufficient accuracy 
in the axial calculations to see error changes for different 
crossflow configurations. There are nine such cases listed in Table 
5. All results have been calculated along a 0 = 33.75° meridian 
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line on the body surface. Since the root mean square error in Cp 
at 0 * 33.75° is very small, for m t choices of panelling and 
inset, detailed plots of singularity strength and Cp distribution 
have not been presented. It was also found that the root mean 
square error In Cp values at other 0 values were essentially 
identical to those presented in Table 2, and that the singularity 
strength distribution solutions were unchanged by 
different 0 meridian line locations of the control points. 

2. 4. 2.1 Singularity type* There are only two choices for 
singularity type In the cross flow calculation, either constant 
doublet or linearly varying doublet. Again, it was found that the 
linearly varying doublet yields better results. Comparing 
configurations 1 and 5, the root mean square error Is reduced from 
0.019056 to 0.005749. Comparisons between configurations 2 and 6, 
as well as 4 and 8, result in the same conclusion. 

2. 4. 2. 2 Panelling scheme. By comparing configurations 1 and 3 , 
it is again found that cosine spacing is better than equal spacing 
for solution accuracy. Reductions -of root mean square error also 
appear in the test pairs configurations 2 and 4 or 6 and 8. 

2. 4. 2. 3 Inset panelling. The Inset panelling technique also 
improves the solution accuracy in the cross flow calculation. This 
technique is recommended for use for blunt body calculations under 
the present singularity panel method. Note that if no inset is used 
for the ellipsoid at angle of attack (not the recommended panelling 
scheme), use of cosine spacing of the cross flow linearly varying 
doublet actually degrades the solution accuracy (see configurations 
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5, 7). This is the only exception in consistency to the conclusions 
made previously for the axial flow analysis. The reason for this 
different trend is not known. 

2. 4. 2. 4 Optimal solution. To find the best panelling arrangement 
to represent the inclined flow past a SR=5 ellipsoidal body, we use 
the optimal case for axial flow, combined with the linearly varying 
doublet with cosine spacing and 1% inset panelling for the cross 
flow. This is configuration 9 in Table 5 which does yield the 
minimum error of all cases tested. Based on these results, use of 
this panelling scheme is strongly recommended for blunt body 
calculations. Figures (8 through 12) display the singularity 
strength and Cp distributions of this optimal panelling arrangement 
at several additional values of a between 5° and 30°. The numerical 
solutions are in excellent agreement with the exact solution 
everywhere on the body surface. 
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CHAPTER 3 


VORTEX LATTICE SURFACE SINGULARITY METHOD 
FOR ARBITRARY BODY OF REVOLUTION 

3.1 Literature Survey 

The surface singularity method is one of the most reliable 
methods available to compute the potential flow past a body of 
arbitrary geometry. By utilizing this method, which is based on 
application of Green's theorem, the original three dimensional flow 
problem is reduced to an integral equation over the given boundary 
surface of the flow region, comprised of the body surface and/or 
wake. The Green's theorem expression represents the disturbance 
velocity potential <j> at any point P in the flow field as the 
integrated effect of source and doublet (or equivalently vorticity) 
singularities on the body surface. The mathematical statement is 

< 1 >(P) = ’ // (7 9 rad ♦ ' ♦ 9 rad 1) • n ds (3.1) 

s 

To reduce this expression to one which contains either source 
only or doublet only, consider a second function, <t>^, which is 
harmonic inside the body, and adding the Green's theorem expression 
for <j>.j evaluated outside the body, one obtains 

( P ) = // Cp- grad (<j> - <j>-| ) ~ (♦ ~ ^ )grad -p-3 • n ds (3.2) 
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Choice of a unique i*-, solution such that <j>, - <p on s results in 

3*1 a4 

a source only method, while fixing <j>, such that -r— = -r£ on the 

I 3n 9n 

boundary yields a doublet singularity method. Also possible are 
various combinations of source and doublet singularities. 

Mathematically, these singularity methods will yield equivalent 
results independent of which representation may be used. However, 
numerical differences do occur when different schemes are used. For 
example, sources are more effective near stagnation points, and 
doublets are more effective in generating and controlling surface 
tangential velocities, since flow near a stagnation point is source- 
like. Intuitively, then, one might expect some benefits to be 
possible from use of the combined source-vortex or source-doublet 
methods. This idea has been investigated and proven correct in 
numerical experiments (Bristow, 1976). 

In summarizing related previous work, discussion of two main 
areas will be given because these two approaches have each had a 
profound influence on the present study. In fact, the computer 
program which has been developed herein to implement the surface 
singularity method is essentially the same as the method due to 
Asfar et al . ( 1 979 ) described below, except for some changes in 
calculation methods and program capability. 

The first approach, due for example to Hess and Smith ‘(1962, 

1967) utilized a source density distribution over the body 

surface. Application of the no-penetration boundary condition on 

the body surface yielded a Fredholm integral equation of the second 

kind for the source density. The profile curves defining the body 

were approximated by a large number of small straight or curved line 
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rectangular elements. The source density was then assumed constant 
over each element. Thus, the contribution of each element to the 
integral expression for the velocity could be obtained by taking 
this constant but unknown source density out of the integral and 
then performing the integration of known geometric quantities over 
the element. Satisfying the integral equation at a control point on 
each element led to a system of linear equations for the source 
densities. Once these densities were evaluated by solving this 
linear system, straightforward substitution in the appropriate 
velocity expressions yielded the velocity distribution and hence the 
pressure distribution. 

The second technique (see, for example, Kandil et al 1974, 
1977; Atta, 1978; Asfar et al.,1979) utilized surface distributions 
of vorticity, in the form of a vortex lattice, and sources, to 
represent the body. The surface source distribution used by Asfar 
et al* (1979) was predetermined in such a way that the source on any 
given element, when acting alone, generated a velocity field which 
cancelled a prescribed fraction of the normal component of the free- 
stream velocity at the control point of that element. The vortex 
lattice was formed by a constant strength quadrilateral vortex 
element, except at the nose or tail where triangular elements were 
used. Description of the formation of the vortex lattice have been 
given by Crigler (1957), Lamar (1971), and Kandil et al*(1974, 

1977). Using the element loop circulations (which automatically 
satisfy the spatial conservation of circulation) as the unknowns, 
and satisfying the no-penetration boundary condition for all 
elements on the body surface, a strongly diagonal influence 
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coefficient matrix of linear equations was generated. Then the 
unknown strengths of circulation for each element were obtained by 
solving this linear system. Since the governing equation is 
elliptic, all aerodynamic parameters were calculated by summing up 
the individual contributions from source singularities and vortex 
lattice circulation elements. Atta (1978) has studied various 
problem formulations in detail, and found the loop circulation 
formulation to be much more convenient than a branch circulation 
formulation. 

Both methods have been found to yield fairly accurate 
solutions, by employing a large number of panels, with a resulting 
high computational cost. According to Asfar et al. (1979), for 
blunt bodies, the second approach of using a combination of source 
and vortex singularities was superior to either the vortex lattice 
only (a special case of the second method, when the source strength 
is zero) or the source panel method; while for slender bodies, the 
vortex lattice method appeared to be superior to the source because 
of the slenderness of the body. This assertion provides another 
justification for using the vortex lattice method in the present 
surface singularity theory. 

3.2 Formulation of the Problem 

The problem formulation for the present method is identical to 
that given in Chapter 1, except that grad F has been replaced 
by n, the outward unit normal vector on the body surface. This 
formulation is more convenient for complicated bodies for which an 
analytical expression for F is unknown. 
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3.2.1 Numerical Scheme 


The body under consideration has been located in a right handed 
cartesian coordinate system with the X axis oriented in opposition 
to the oncoming flow (see Sketch 3.1). The undisturbed free-stream 
is, in general, oriented at an angle of attack, a, with respect to 
the X axis, in the XZ plane. The body has been described by a 
number of discrete points on it's surface. In the present work, 
where the surface of the body was described by an analytic function, 
these body points have been generated by the program itself. For 
more complicated geometry, the program user must define these body 
points manually. 

The surface points have been organized in "rows" and "columns" 
as indicated in Sketch 3.1. The lattice consists of short straight 
lines connecting the body points, thereby forming quadrilateral 
elements, with the exception of the nose and tail, where these lines 
converge so that the elements are triangular. A system of elements 
which are chosen to be quadrilaterals of aspect ratio of nearly 
unity has been found to yield the best results. The velocity fields 
generated from the above scheme satisfy equations (1.4 and 1.6) 
regardless of the circulation strengths. It only remains, 
therefore, to choose their strengths so that equation (1.5) will be 
satisfied. Equation (1.5) can also be manipulated into a set of 
algebraic equations, by substituting the following expressions into 
the equation: 


7 *i " j S 1 ( V/ + q u r / + r u r j l?) 


(3.3) 
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"j 58 l + B.j + 


(3.4) 


= ’ U -COSat - U^slna t (3.5) 

where the subscript i denotes the associated quantity evaluated at 
the control point "1% T. Is the unknown strength of circulation of 
element "j"; P^., q^, r^ are the Influence coefficients, which 
give the velocities In the X, Y and Z directions at the 1 th control 
point due to a unit loop circulation around the j th element, 
and A f , B. , are the direction cosines of the surface unit normal 

at the control point "i". Expressions for the influence 
coefficients will be given subsequently. 


Substitution of equations (3.3 to 3.5) into equation (1.5), 
yields a set of equations for the unknown loop circulations as: 



Sketch 3.1 Coordinate scheme for the problem 


A i j!;, Vj + B t ^ Vj + c i j f, r ij r j 


= Ujrosa A. + u^slno for 1 = 1,2, ... N (3.61 
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The sequence of numerical techniques which are involved in 
solving this linear system include decisions as to how to 
efficiently store the large influence coefficient matrix, selection 
of a linear equation solver, and reduction of the size of the linear 
system due to symmetric properties, as has been well documented by 
Kandil (1974) and Asfar et al. (1979). The following two sections 
are devoted to development of the explicit expressions for the 
constant coefficients which appear in equation (3.6). 

(x 2 ,y 2 , 

X 

Z 

(a) A closed loop vortex element 


(*i >y-| »^i ) 

(b) A finite length vortex filament 
Sketch 3.2 Nomenclature for calculating the induced velocity 
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3.2.2 Calculation of Influence Coefficients o. . ,_q. r. 

T J 1 j T J 

The total induced velocity at the i th control point due to a 
discrete vortex element "j", with the constant strength of 
circulation along the sides (see Sketch 3. 3. a), can be expressed in 
the following symbolic form: 

^induced “ P i j r j 1 + q ij r j J + r ij r j k 


4 


= z 
4=1 


(u 1t r j T + v u r jJ + w U r j (?1 


(3.7) 


The induced velocity due to a single finite vortex segment (see 
Sketch 3.2.b) is calculated by using the Biot-Savart law (see 
Karamcheti, 1966). 


r • e 2 r . 

tl = writ, s1ne d 9 e, * Tif- ,C0S9 1 - =»s8,)e (3.8) 

PI D 


„ I ? 1 x \\ 


r • r * > 

where r_ * — , cose, = 1 , cose, = V/ 2 (3.9) 


Vl 


r o r 2 


and e y , the direction of the induced velocity, is given by the unit 
vector 





(3.10) 


substituting equations (3.9 and 3.10) into equation (3.8), a vector 
expression for this induced velocity can be obtained as follows: 


fr = J 
ll Tit 
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(3.11) 



where r Q = (x 2 - x^t + (y 2 - + (z 2 - )it 

r 1 * (x c - x 1 )i + (y c - y-i )j + (z Q - z ] )£ 


r 2 = (x c - Xg)t + (y c - y 2 )j + (z c - z 2 )ic (3.12) 

^2’ ^l3* are ca l cu ^ atecl * n the same way- By summing the 
contribution of each segment together, P^, q^j, and r^ are given 

T 1 ? (y c ' V (z c - Vl> - ,y c - Vl Hz c - z *> 

« * jt»i *Ti 


(W t2 - w l3 ) 


q 1j 


w *=i 


(x c ' Vl )(z c - z *> • (x c - 
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(W , 2 . W, 3 ) 


r 1 j = W Z 


1 t (x c ' x t )(y c -W - (x c - 


VI 


*= 1 


'*1 


)(y c - y.) , 

(W *2 " W *3 ) 
(3.13) 


where 

“j-l = C(y c - V (z c • Vi 1 - (y c - y i + l )(z c - z i )]2 


+ [(x c - x » ),z c ■ W - (x c - x «H )(z c - z « )]2 
+ C(X C - x t My c - y Itl ) - (x c - x 1+1 )(y c - y ,)] 2 (3.14) 

u _ . ( Vi - x « 1(x c - V - (y tn - y t »y c - y t> * (z t+i - z t )(z c - h> 

V< x c - x t ) 2 * (y c - y/ + (I c - 3 / 

( Vl ~ x t l(x c - W * (y ui - V ty c - W * (z x+l ~ z * )fz c - W 

V (x c ■ x t+i* + ^ + ^ 2 c “ z j+l ' 2 51 
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These expressions have been implemented in a FORTRAN computer 
program to compute the velocities due to a closed loop discrete 
vortex element. Note that the index of t+l is set equal to 1 
when i is 4. 


3.2.3 Unit Normal Vector Calculation and Control Point Selection 
The unit normal vector and the control point coordinates can be 
obtained through algebraic combination of the four given element 
corner points (see Sketch 3.3), as is now described. 



Sketch 3.3 Nomenclature for calculating the unit normal vector 
In order to simplify the geometry calculation, it will be 
assumed that the element is planar; that is, the four points 
defining a panel element all lie in the same plane. This assumption 
has been satisfied by the axi symmetric bodies studied in the present 
work. 


3. 2. 3.1 Unit normal vector calculation. From equation (3.4), the 
unit normal vector Is defined as n^ • f + B. J it, where A^, 
B.j, and are constant coefficients which are determined by the 


- r l x r 2 

element geometry, n. can be represented by - — , where 

' \\ * r 2 \ 
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?1 = (x 2 - x 4 ) f + (y 2 - y 4 ) j + (z 2 - z 4 ) It 
? 2 = (x 3 - x t ) t + (y 3 - j + (z 3 - z^ It 
Therefore, 

x r 2 = C(y 2 - y 4 )(z 3 - z^ - (y 3 - y 1 )(z 2 - z 4 
+ C(x 3 - x 1 )(z 2 - z 4 ) - (x 2 - x 4 )(z 3 - z x )] j 

+ C(x 2 - x 4 )(y 3 - y x ) - (x 3 - x x )(y 2 - y 4 )] k* 

= A' t + B' j + C' t 
By definition, then 

A. = £ 

"j r— — 

V (A * ) 2 + (B ‘ ) 2 + (C ' ) 2 

B i B> 

V(A ') 2 + (B * ) 2 + (C * ) 2 

c = C' 

where (B') 2 + (C ) 2 

A ' * (y 2 - y 4 )(z 3 - z x ) - (y 3 - y^iZg - z 4 ) 

B' = (x 3 - x^izg - z 4 ) - (x 2 - x 4 )(z 3 - Zj) 




( 3 . 15 ) 


T 


( 3 . 16 ) 


( 3 . 17 ) 


C ~ (x 2 - x 4 ) (y 3 - y-j) - (x 3 - x^)(y 2 -y 4 ) 


( 3 . 18 ) 
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The choice of location of the element control point, where the 
no-penetration boundary condition is satisfied, is crucial due to 
the singular nature of the discrete vortex filaments. Following 
Hess and Smith (1962), Asfar et al. (1979) concluded that the best 
results were obtained when the control points were located so that 
the normal component of self-induced velocity due to the closed loop 
vortex filament at it's own control point is a minimum (known as the 
"null point" according to Hess and Smith). These minimum self- 
induced velocity points were located through a tedious search 
technique. The search scheme is costly and is restricted to 
problems of fixed geometry. For panels with edge separation (for 
example in the wing problem), this has to be repeated for the same 
panel at each iteration. Hence it is not a practical scheme for the 
wing-body combination problem. Therefore, in the present study, 
this scheme has not been used. Instead, the centroid of the panel 
area has been chosen as the control point location, since for the 
axi symmetric bodies studied, all panels are planar. The advantage 
of this choice of location is the ease of calculation of the 
centroid, provided there is not too much difference for the control 
point location between these two modes. Numerical experiments have 
proven that this difference is small for both rectangular and 
triangular elements. The details of calculating the centroid and 
the panel area for a typical quadrilateral element is given in 
Appendix 8. Also shown in Appendix 8 are the results of the 
numerical experiments comparing self-induced velocities for various 
control point locations. 
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3.2.4 Pressure Coefficient Calculation 


Once the linear system given by equation (3.6) has been solved 
for the unknown circulations, the pressure coefficient at each point 
on the surface of the body is calculated by using Bernoulli's equa- 
tion. 


c n = 1 - -X (3.19) 

P 7? 

00 


Here, Y t is the total tangential velocity of the flow at the control 
point, and is evaluated, as described in Kandil et al. (1977), as the 
sum of 3 parts; first, the part due to the influence of all the 
vortex loop elements in the lattice at the control point, second, 
the contribution of the free stream velocity, and third, the 
velocity jump across the vortex sheet. The third term is obtained 
by estimating an equivalent continuous vortex sheet from the lumped 
vortex filaments. This contribution is to account for the self- 
induced tangential velocity due to the local strength of the vortex 
sheet as explained by Kandil (1981). If this term is not taken 
into account, errorous pressure coefficients are found. For panels 
with continuous vorticity, this term automatically exists in the 
induced velocity expression (Kandil et al., 1980). Kandil et al. 
(1977) have given a detailed description for a flat, rectangular 
element. Mote that there is a typographical error in the definition 
sketch in that paper. Extension of this method to a general, 
quadrilateral vortex element has been implemented by Atta (1978) for 


*Private discussion with Dr. Kandil (1981). 


55 




a control point located at the average of four corner points. These 
two techniques for calculating the tangential velocity have been 
summarized in Appendix C. In the present studies, a modified method 
has been developed to calculate the tangential velocity due to the 
local quadrilateral vortex element. The idea is based on the 
equation due to Kandil et al . ( 1 977) then extended to a quadril ateral 
element. The reason for developing the new expression instead of 
using the formula due to Atta (1978) is the lesser generality of the 
latter expression. It is easy to find that the formula due to Atta 
(1978) will give the same magnitude of tangential velocity for any 
point on the element, which is not consistent with the mathematical 
model. The modified method has the ability to account for the 
deviation of the tangential velocity between different points. For 
some special cases, as shown in the next section, the modified 
method gives solutions which are identical to those obtained using 
the expression due to Atta (1978). The comparison of these two 
schemes for either axi symmetric flow or inclined flow will also be 
shown in the next section. 

For a typical planar quadrilateral vortex element with two 
parallel sides, having the other two sides of equal length, as shown 
in Sketch 3.2, the tangential velocity of a point located on the 
center line is given as follows: 


1 1 


T vx " TIT 


2a, 


Cr 'l T d 2 + r '3 a d l-* 
P P 


£ i 

T vy = [r 4 cos8 / + r 2 cos0 T* 


2a 


(3.20) 
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d. d.tane d. d.tane 

Where ■ r t ♦ r-,sine - —— - lylne - — 


d 2 d 2 tan 0 d 2 d 2 tan 0 

r ' 3 * r 3 * r 2 sine q ~q~ - V 1 " 9 q ~q~ 



Sketch 3.2 Nomenclature of calculating the induced 
tangential velocity 


These expressions have been developed from those of Kandil et al. 
(1977) (summarized in Appendix C), and are similar to those of 
Kandil et al. (1982). Modifications have been made to account for 
the nonparallel orientations of r 2 and r 4 , as well as to account for 
the unequal lengths of filaments and r^. The ratios z^/z and 
43 /Zp in the r yx expression are to account for these unequal fila- 
ment lengths. The filaments r 2 and r 4 have been resolved into com- 
ponents along the x and y axes, and their x components used in 
computing T yy . The signed y components of r 2 and r 4 have been used 
to replace and by rj and as given above. Note that 
these r 2 , r 4 y components have been lumped with and r 3 using 
linear interpolation, to allow movement of the control point away ^ 
from the panel center. 


In case the field point is not located on the body surface, but 
is very close to the vortex segment, there is no tangential velocity 
component in the total velocity expression. However, the induced 
velocity due to the vortex segment must be modified by a viscous 
core relation in the form of an exponential multiplier, which causes 
the induced velocity to approach zero as the vortex is approached. 
The reason for this is simply because the law of Biot and Savart 
which is based on an inviscid flow, is singular. However, no arti- 
ficial viscosity coefficient has been added in the present 
technique. An accurate accounting of the viscous effects must be 
obtained from boundary layer solutions, which is beyond the scope of 
the present study. Instead the velocity and pressure coefficients 
only at points on the body surface or at far field points have been 
obtained in the present work. 

3.3 Accuracy of Method 

Since the surface panel method has been developed for solving 
the potential flow past an arbitrary body, therefore, the technique 
of using the vortex ring element as an option to solve for the 
problem of axial flow past an axi symmetric body of revolution has 
not been considered here. Note that a very accurate C p distribution 
could be obtained using this technique by employing a large number 
of ring elements. The test configurations investigated in the 
present study still remain in the axi symmetric body category. 

However, the flexibility for modelling a realistic fuselage 
(symmetric with the y=0 plane) has been considered in the computer 
program. Because of the assumed symmetric property, only one half 
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of the control points on the body surface are stored, but the 
influence coefficients of the induced velocity at these points are 
evaluated in the usual way. For each one of these control points, 
the influence coefficient due to the image of the stored element is 
evaluated by simply reflecting the coordfnates of the corners of 
this element. The contributions are then added to those of the 
stored element. This reduces the number of unknowns by one half. 

In studying the performance of the method for the test cases 
discussed below, a serious problem of the need for a large number of 
elements to adequately represent the given body has been found to 
exist. This situation is undesirable in terms of cost. However, 
the reliability of the present scheme may be recognized by studying 
the solutions presented below. 

In this chapter, as in Chapter 2, the test geometry which has 
been considered and compared with the corresponding exact solution 
is a SR=5 ellipsoid. The tangential velocity jump has been 
calculated by the three different methods listed above: Kandil et 

al . (1977), Atta (1978), and the modified method, equation (3.20), 
and results have been displayed and discussed in sequence. 

3.3.1 Axi symmetric Flow 

For calculation of the axi symmetric flow past an ellipsoid of 
SR=5, a panelling scheme using cosine axial distribution of panels 
has been applied, to Improve the solution accuracy. This is 
consistent with the criteria mentioned in section 3.2.1, where it 
was Indicated that the best results would be obtained if the vortex 
panel elements were chosen to be of aspect ratio of nearly unity. 

Of course, the cosine spacing axial distribution is not the only way 


to represent the body. The optimal panelling technique Is expected 
to vary as the body geometry is varied. Figure 13 shows that the C p 
distribution for the SR=fT, ellipsoid along a meridian line obtained 
by using 8*30 panels (8 panels in the 0 direction by 30 panels 
axially), with a cosine spacing axial distribution, is much more 
accurate than that obtained using the same number of panels but with 
equal spacing in the axial direction. 

In order to discover the differences among the three tangential 
velocity calculation methods, Table 6 has been prepared. The test 
runs have all been based on the same panelling scheme as above of 
8*30 panels, with a cosine axial distribution, and the control 
points have been located at the average of the four corner 
coordinates of each panel. In fact, location of the control points 
at the average of the four corner points is not recommended in 
general calculations. Another comparison, based on the recommended 
control point location of the centroid of area of each panel, will 
be given below. Note that the pressure distribution obtained by 
using the method of Atta (1978) is exactly the same as that which is 
obtained by using the modified method (equation 3.20), provided that 
the control points are selected to be the average of the four corner 
coordinates of the panels. Atta's application of Kandil's method 
and the present modified method both work better than Kandil's 
original formulation does, for the non-rectangular elements utilized 
in this test case. Especially at the nose and tail region, where 
the element is less rectangular, the solution accuracy has been 
improved by accounting for the length effect and orientation of the 
vortex segments. The same conclusion is reached 


60 


in the second comparison, as shown in Table 7, where the control 
points were selected to be the centroid of area of each panel. 
Another way to improve the solution accuracy, by moving the control 
points of all triangular elements from the centroid of area to the 
3/4 chord length away from the apex of the element, was attempted, 
listed as "Scheme 2" in Table 7. No significant improvement was 
achieved through this effort. Numerically, Atta's method appears to 
be slightly superior to the modified method. However, the observed 
differences between these two methods are judged to be negligible. 
Therefore, the modified method (equation 3.20) is recommended in 
future applications, since it correctly accounts for variation in 
control point location. It is also noted that the solutions 
obtained from the vortex lattice surface panel method are less 
accurate than those obtained from the axial singularity method. The 
reason is because the present calculation doesn't involve enough 
panels to adequately represent the body, and due to the poor 
accuracy obtainable with the discrete vortex singularities used to 
represent the continuous surface vortex sheets. 

3.3.2 Inclined Flow at 20° Angle of Attack 

Procedures identical to those used in the axi symmetric flow 
analysis above were repeated for the flow at a 20° angle of attack 
for the SR=5 ellipsoidal body. Table 8 shows the C p distribution 
along two meridian lines, one located 7.5° away from top meridian 
line, with the other located 82.5° away from top meridian line, 
using 12*20 panels, cosine spacing axial panel distribution, and 
with the control points located at the average of the four corner 
coordinates of each panel. Numerical solutions obtained from both 


Atta's and the modified method match the exact solution fairly well 
along the 0 = 82.5° meridian line, but are less accurate near the 
nose along the 0 = 7.5° meridian line. This may be explained in 
that the assumption of a linearly distributed vortex sheet strength 
used in the tangential velocity calculation is no longer adequate at 
the leeward side of the body near the nose when a large angle of 
attack is imposed. 

Tables 9 and 10 display the C p solutions obtained by moving the 
control points to the centroid of area of the panels and keeping the 
rest of the conditions the same as for Table 8. Through comparison 
of these results, it is found that taking the centroid of the panel 
area as the control point location yields a more accurate C p 
solution. This reinforces Asfar et al.'s assertion about choosing 
the minimum induced velocity point as the control point location to 
obtain the best solution. 

Figure 14 displays the C p distribution along four meridian 
lines in the first quadrant, compared with the exact solution, for 
the SR=5, a = 20° ellipsoidal body, using 12*20 panels, cosine 
spacing axial distribution, and the modified tangential velocity 
calculation method. Figure 15 displays the C p distribution along 
six circumferential lines at different x/L locations and compares it 
with the exact solution for the same body using the same method. In 
general, the calculated solution matches the exact solution well, 
except at the leeward side of the body near the nose. Note that 
the 0 coordinate used in the surface singularity program is shifted 
180°_compared to that used in Chapter 2 for the axial singularity 
method. 
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CHAPTER 4 


RESULTS AND COMPARISONS FOR GENERAL 
AXISYMMETRIC CONFIGURATIONS 

4.1 Introduction 

To further investigate the accuracy and the capability of the 
two computer programs described in Chapters 2 and 3, six additional 
configurations have been investigated in a series of numerical 
experiments. The main area of concern for the axial singularity pro- 
gram has been to investigate the feasability of obtaining 
approximate C p distributions for a series of axi symmetric bodies 
having discontinuous surface curvature, such as a cylindrical body 
with an ellipsoidal nose and tail. Previously, Zedan and Dalton 
(1980) have claimed that the axial singularity method cannot 
accurately model bodies with either discontinuous surface slope or 
surface curvature. It has been found in the present work that the 
axial singularity method appears to have the capability of modelling 
the second type of configuration having discontinuous changes in 
surface curvature, provided a careful numerical scheme is applied, 
as described below. Solutions for these types of configurations 
using the axial singularity method have been obtained and compared 
in this chapter with those obtained from the surface source program 
of Hess and Smith (1967) and available experimental data. A brief 
discussion about the panelling techniques has also been given 
separately for each test case. The vortex lattice surface 
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singularity program described in Chapter 3 has been run for each 
configuration, and compared with all other theoretical and 
experimental results. An explanation has been given in the 
discussion of the last configuration (30® cone-cylinder), to clear 
up the faulty conclusion which appeared in Asfar et al. (1979), where 
it was claimed that the surface vortex lattice method was superior 
to the surface source method. Generally, the comparisons shown in 
this chapter demonstrate that good Cp distributions can be obtained 
from either the axial singularity or the vortex lattice surface 
methods. Despite the fact that the surface velocity is extremely 
sensitive to the details of the surface geometry, the above two 
methods obviously have the capability to account for a wide range of 
surface geometry accurately. This means that the solution accuracy 
is largely dependent upon the input surface geometry and panelling, 
but not the methods themselves. This chapter is intended to serve 
as an indication of the capabilities of the two methods, and is 
therefore not intended to be all-inclusive. In particular, no 
results have been obtained for general three-dimensional bodies. 

This is largely due to a lack of effort to develop or utilize a 
general three dimensional geometry input routine such as GEMPAK (see 
Stack, 1982). 


4.2 Test Configurations 

4.2.1 Modified Ellipsoid with Zero Curvature at the Maximum 
Thickness 

This body is similar to a configuration which was described by 
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Hess and Smith (1967) as a two-dimensional body of SR=5, whose 
profile curve was given by 


y = + 0.2^ 1 - x 3/2 (4.1) 

Despite the strong similarity in the shapes of the same slenderness 
ratio 5 ellipsoidal body of revolution and this configuration, their 
surface pressure distributions are, of course, quite different. 

That for the ellipsoidal body is very flat and entirely concave 
down, while that for the modified ellipsoidal body obtained from 
equation (4.1) by replacing y by r, has a sharp peak at the shoulder 
and a wide region where it is concave up. 

The reason that this configuration has been chosen as a test 
case in the present study is that the geometry contains an infinite 
radius of curvature at the shoulder, and as a result, it can be 
viewed as a transitional model from the ellipsoid to a‘ cylindrical 
body witn an ellipsoidal nose and tail. 

Figure 16 (a) shows the source strength distribution along the 
body axis obtained by using the constant and linearly varying axial 
singularity schemes described in. Chapter 2. Both cases used 20 
cosine type discrete panels to represent the given body. The 
difference in the strength distribution for the linear source 
singularity method between the ellipsoid and this modified ellipsoid 
configuration is that the former one (Figure 2) has a continuous 
strength at the shoulder and the latter one has a jump in strength 
there. Apparently, the zero curvature at the shoulder is 
responsible for this jump. One must account for this discontinuity 
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in source strength when modeling those classes of bodies with 
discontinuous curvature using the linear source singularity. In 
building up the influence coefficient matrix, the continuity 
equation at the shoulder must be replaced by a constraint equation 
which will allow a strength jump but having zero net strength 
there. Therefore, the anti -symmetric property in axial singularity 
strength distribution is still maintained without any change due to 
the effect of discontinuous curvature at the shoulder. Figure 16 
(b) shows the C p distribution along a meridian line at a = 0°. Both 
constant and linearly varying source appear to give the same level 
of solution accuracy. The C p curve also has a cusped peak at the 
shoulder which is the characteristic C p property for this geometry. 

In order to confirm the accuracy of the solution obtained from 
the axial singularity program, the same configuration has been 
analyzed using the vortex lattice surface singularity program. The 
panelling scheme undertaken herein was an 8*30 cosine type axial 
distribution. Figure 17 displays the comparison of C p distributions 
between the two different methods. The equivalence of the two 
solutions shown does provide evidence that the two solutions are 
each correct. No experimental data for this type of body is known. 

4.2.2 Ogival Body of Revolution 

The surface of this configuration is composed by rotating a 

finite circular arc around the body axis. Such geometries have a 

pointed nose and tail. The body surface is governed by an analytic 

function where the slope and curvature are continuous. At the 

region very close to the pointed nose and tail, potential flow 

theory will not predict an accurate C solution. However, it is 
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believed that a satisfactory C p distribution can be obtained by 
using the present two programs over most of the body surface. 

Figure 18 (a), (b) display the source strength distributions 
along the body axis obtained using constant and linearly varying 
axial singularity schemes for SR=4 ogival body of revolution. The 
panelling technique undertaken for both computer runs is 20 discrete 
cosine spaced panels, without inset, to represent the given body. 

The selection of zero Inset distance has been made due to the fact 
that the radius of curvature of a pointed nose is zero. This is 
consistent with earlier recommendations of Zedan and Dalton 
(1980). Note that the source strength at the nose and tail is found 
to be zero, and the strength distribution is strictly anti- 
symmetric. The appearance of these two properties is helpful in 
judging the solution quality. Figure 18 (c) displays the C p 
distribution along a meridian line. The constant and linearly 
varying source schemes give nearly identical C p values for the 
entire region except at the nose and tail where the pressure has a 
sudden drop if the former scheme is used. The explanation for this 
incorrect pressure drop is that the constant strength axial 
singularity distribution is Inadequate to describe flow at the nose 
region of an ogival body. However, the solution accuracies away 
from the nose and tail are strongly confirmed in this comparison. 

As was done in analyzing the first configuration, Figure 19 has 
been prepared to show the comparison of C p solutions between the 
axial singularity and the vortex lattice surface panel programs. 
Again, a panelling scheme using 8*30, cosine axial spaced panels has 
been applied to model the present configuration in the vortex 
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lattice surface program. The discrepancy appearing In this figure 
is believed to be due to error in the surface vortex lattice 
method." An equal spacing panelling scheme reduces the numerical 
error (see Figure 19), but perhaps addition of some source 
singularities neqr the stagnation points may be necessary to obtain 
more accurate results. 

4.2.3 Cylindrical Body with Modified Ellipsoidal Nose 

This configuration consists of a cylindrical middle section 
with the radius equal to the semi -mi nor axis length of the modified 
ellipsoidal nose and tail. Therefore, the surface slope is 
continuous everywhere and surface curvature goes smoothly to zero on 
the middle cylindrical surface. This geometry is a special example 
of a smooth body of revolution. Slender body theory can well 
predict the C p distribution in the middle part of such bodies, but 
not for the nose and tail. The present study does not impose the 
slender body assumptions in formulation. As for the first 
configuration, a strength jump at the junction where the radius of 
curvature goes to infinity will be allowed when the linearly varying 
source distribution is used. 

Figure 20 (a), (b) show the axial source strength distributions 
along the body axis obtained using constant and linearly varying 
source distributions. Both cases use 10 panels for the nose and 
tail, and one panel for the cylindrical part, where the cross 
sectional area is constant. In the case of using linear type 
singularity distributions to model the given body, the two source 
strength continuity equations at the junctions have been replaced by 
two extra no-penetration boundary equations. In other words, two 
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extra control points have been located on the body surface. 

Numerical experiments show that the Cp distribution in this 
configuration is not sensitive to the location of these two extra 
control points. However, some adjustment of control point location 
is needed for similar bodies having a long cylindrical middle 
section. More discussion about this adjustment will be given in 
discussion of the fifth configuration. 

The Cp distributions from these two test runs are nearly 
identical, as shown in Figure 20 (c). Although there is no 
experimental data available, this configuration has been modeled 
using the vortex lattice surface panel program. Panelling used is 
8*30 panels with equal axial spacing. Through the consistency of 
the Cp solutions for the above three theoretical cases, one may 
conclude that both the axial singularity method and the surface 
vortex lattice method are adequate for this configuration. 

4.2.4 Cylindrical Body with Ogival Nose 

This configuration was obtained by replacing the SR=5 modified 
ellipsoidal nose and tail with a SR=4 ogival nose and tail. The 
detailed geometry is shown in Figure 21. A similar body with a boat 
tail, as shown in the same figure, was utilized in a wind tunnel 
experimental test by Fox (1971) to obtain Cp distributions at a Mach 
number equal to 0.403. Over the front half of the body, the 
pressure distribution will not be strongly influenced by the 
different tail shapes. Therefore, the same symmetrical 
computational model has been used for both the axial singularity 
program and the vortex lattice surface panel program. Good 
agreement with the experimental Cp data is observed in Figure 21 for 
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all methods. The solution obtained from the axial singularity 
program matches the experimental data even better than that from the 
vortex lattice surface panel program. Of course, the small number 
of vortex lattice panels used in modelling this configuration is 
likely responsible for this discrepancy. 

4.2.5 Cylindrical Body with Ellipsoidal Nose 

This configuration is composed of a SR=2 blunt nose and a long 
cylindrical body. The available experimental Cp data were obtained 
by Campbell and Lewis (1955). The panelling technique used in the 
previous two configurations has been found to be inadequate for the 
present configuration. The reasons are, first, the nose is too 
blunt, and second, the length of cylinder section is too long. It 
has been found that the Cp solution for this configuration is very 
sensitive to the panelling and control point location, a result 
which was not seen for the previous two cases. In general, it has 
been found that the constant source singularity scheme cannot handle 
geometries similar to the present configuration. Based on 
accumulated experience, a new panelling scheme has been discovered 
and found to be adequate to yield a set of reasonable C p values, for 
the linear source scheme. Again, the discontinuous strength at the 
juncture is allowed. Also, three panels are used to represent the 
cylindrical section. Since a sudden localized sink distribution was 
expected to be necessary to accurately model the flow around the 
shoulder over the junction between the nose and cylinder, the 
additional two panels on the cylinder were chosen to be relatively 
small. Because the two source strength continuity equations at the 
junctions do not exist, two extra control points were needed to 
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perform the calculations. Numerical experiments have shown that 
these two extra control points, combined with the original three 
control points, corresponding to the three panels in the cylindrical 
section, should be equally distributed on the cylinder surface to 
yield the best results. The detailed control point locations on the 
cylindrical section have been shown along with the resulting C p 
solution in Figure 22. It has been found that the size of the two 
small panels is not critical, between one percent and four percent 
of the cylinder length, in obtaining the same accuracy of Cp 
distribution, even though the singularity strength distribution is 
changed somewhat. Figure 22 shows the Cp distribution along a 
meridian line when the flow is at zero angle of attack. Comparing 
the solution obtained from the axial singularity program to 
experimental data, a good agreement is noted over the front portions 
of the computational model. The small discrepancies on the 
cylindrical surface are mainly due to differences in geometry 
between the wind tunnel model and the computational model. This 
assertion has been verified by calculating the Cp values for a 
configuration having a longer cylindrical middle section, with the 
same size nose and tail. Pressures for such a geometry were closer 
to the data. Figure 22 also shows the Cp solution obtained from the 
vortex lattice surface panel program. The panel scheme used herein 
is 10*32 panels with a cosine type axial distribution. The solution 
accuracy for the vortex lattice surface panel method is seen to be 
not as good as that of the axial singularity method. 

Results are presented in Figure 23 for the same configuration 
as in Figure 22, but at an angle of attack of a = 6.08°, along the 
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top meridian line of the body. Experimental results have been shown 
as a solid curve because no data values were given by Campbell and 
Lewis (1955) when a t 0°. Panelling, inset and control point 
location for the source flow are the same as was used in Figure 22 
at a * 0°. However, the linear crossflow doublet calculations have 
been made assuming continuous strengths at the junctions between 
nose or tail and cylinder. Also, control point location on the 
cylinder has been chosen to be of the cosine type. Agreement 
between experiment and the axial singularity method is not as good 
as was observed for o = 0°; the reasons for the observed 
discrepancies are not known. Also shown are results at e = 9° for 
the vortex lattice surface method, using 10*30 panels; these results 
are less accurate than the axial singularity results. 

4.2.6 30° Cone Cylinder 

This configuration first was studied by Johnson (1963). He 
compared pressure distributions computed by the surface source panel 
method (Douglas Neumann Program) of Hess and Smith (1967) with 
experimental data which he obtained at zero angle of attack and at 
plus and minus 20° angle of attack. The cylindrical after-body of 
the wind tunnel model had a length equal to 2.9 times its radius. 

The geometrical model utilized in the surface source panel program 
was semi-infinite. This difference between the geometries of the 
experiment and theory led to discrepancies between calculated and 
experimental pressures on the after portions of the cylindrical 
body, but its effect was found to be negligible over the nose 
regions. More recently, the same configuration has been studied by 
Asfar et al . (1979) to test the solution accuracy for their vortex 
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lattice surface panel method. The computational geometry in their 
work was exactly the same as the wind tunnel model of Johnson. 
Therefore, the pressure distributions obtained in Asfar's work 
matched the experimental data much better than those obtained from 
the surface source panel program of Hess and Smith over the after 
portions of the cylindrical body. Asfar was apparently unaware of 
these differences In geometry, and mistakenly concluded that the 
surface vortex lattice method was superior to the surface source 
method. In order to show that this conclusion is incorrect, two 
test runs have been prepared to illustrate the effect of after-body 
length of the 30° cone-cylinder geometry. As shown in Figure 24, 
the Cp distribution of the first geometry, having a cylindrical 
after-body with a length equal to 2.9 times its radius, matches the 
experimental data well. When the cylindircal after-body length was 
extended to be 8.5 times the cylinder radius, then the Cp 
distribution of this second geometry matches the solution obtained 
from the surface source panel program. Note that all four pressure 
curves agree well along the cone surface. In the after-region, the 
predicted pressure values could be strongly influenced by the sharp 
corner. 

Figures 25 through 27 compare the calculated and the experi- 
mental pressure distributions along three different meridian 
lines (0 = 45°, 0 = 90°, 0 = 135°) at 20° angle of attack for this 
configuration. At this angle of attack, agreement is good over the 
nose portion of the body for all cases, and less satisfactory over 
the after-body. The biggest discrepancy occurs at the leeward side 
of the after-body, as shown in Figure 25, where the flow is expected 
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to be separated for this angle of attack. 

Figure 28 compares calculated and experimental circumferential 
pressure distributions at two axial locations on the pointed 30° 
cone-cylinder at 20° angle of attack. Again, the effect of 
separation on the leeward side of the cone can be seen, particularly 
for the after station. Panelling used for this configuration is 
10*32 panels, equally distributed along the x axis. Of course, the 
solution may be improved by adequately refining the panels. It is 
worthwhile to note that no success was had in attempts to obtain the 
Cp distribution using the axial singularity program for this 
geometry. 
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CHAPTER 5 


CONCLUSIONS 

Two different types of singularity methods have been considered 
in solving for the potential flow due to a uniform stream, either 
axi symmetric or at an angle of attack, past a non-lifting body. The 
theoretical formulations have been described in separate chapters of 
this thesis. Two FORTRAN computer programs have been generated and 
tested to implement these theoretical models. Solutions obtained 
from these two programs for certain classes of axi symmetric bodies 
have been compared with exact solutions or available experimental 
data. Numerical experiments have shown that both methods give 
accurate Cp distributions for general axi symmetric bodies of 
revolution, except that the axial singularity method cannot handle 
bodies with discontinuous surface slopes. For axisymmetric bodies 
which are sufficiently smooth, the axial singularity method is more 
accurate than the surface vortex lattice method. Some important 
features about these two methods will be summarized below. 

The axial singularity method is considered to be a method which 
is quite inexpensive, since little computer memory is required, but 
it has more restrictions on body geometry compared to surface 
singularity approaches of solving for the potential flow past an 
axisymmetric body. In the course of analyzing an aircraft 
configuration, the axial singularity method will not be adequate for 
modelling a general fuselage. However, at the preliminary design 
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stage to estimate the interactions between a fuselage of simple 
geometry and the non-planar wing, this method is far more applicable 
than any other numerical approach. By employing higher-order 
singularity distributions, great improvements in the solution 
accuracy and geometry flexibility have been found. Specifically, a 
higher-order axial singularity distribution model has been 
formulated for either axial or cross flow, where constant and linear 
orders in strength variation have been utilized to yield a numerical 
solution. It has been found that panel and control point location 
for the higher-order singularity method can be very important for 
obtaining a good quality solution. Inadequate choice of panelling 
for this technique will yield solutions which are incorrect, as was 
found, for example, using the continuous linearly varying source 
distribution applied to model the flow past a cylindrical body with 
either ogival or ellipsoidal noses. Judgement of the adequacy of 
the solution can be aided by visualizing the variation of source 
strength distribution along the body axis. For ellipsoidal bodies, 
a reliable and accurate panelling technique has been described, 
where the optimal Cp distribution for the axi symmetric flow past a 
general ellipsoid has been found by using cosine spaced, linearly 
varying source distributions with the suggested inset distance 
governed by equation (2.69). The solution accuracy so obtained has a 
root mean square error in pressure coefficient on the order 

_5 

of 10 . The same panelling and inset has been found to yield 

comparable solution accuracies for ellipsoids of slenderness ratios 
between 2 and 10 at angles of attack up to 30° degrees. Therefore, 
the inset distance plays a key role in obtaining the optimal Cp 
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solution for blunt nosed bodies. This optimal inset appears to 
depend most strongly upon the local radius of curvature of the body 
nose. A more limited effect of the size of the singularity element 
nearest the nose mey also be noted. This discovery can be viewed as 
being contradictory to assertions made by Zedan and Dalton (1980), 
where it was recommended that 1% inset distance be used for all 
ellipsoidal body calculations. Another important feature in the 
present stu<ly has been the attempt to obtain solutions for 
axi symmetric bodies with discontinuous surface curvature using the 
axial singularity method. Despite the fact that some reasonably 
good Cp distributions have been obtained from such geometries by 
using a specialized singularity panel size variation consisting of 
cosine panel distributions on the nose and tail, plus two very small 
panels at the ends of the cylinder, no general panelling rule 
applicable to all configurations has been developed at the present 
time. 

Accuracy of the vortex lattice surface panel method implemented 
in the present study has been investigated using various 
axi symmetric configurations. The Cp distributions have been 
compared with exact solutions, experimental data, and results of the 
surface singularity theory. It has been found that the solution 
accuracy is largely dependent upon the panel aspect ratio and the 
control point location. Best results will be obtained by chosing 
the panels to have aspect ratios near unity and taking the minimum 
induced velocity point inside the panel to be the control point. 

The panel centroid has been found to be an acceptable control point 
location. The local tangential velocity calculation method 
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devleoped in previous works by Kandil et al . (1977); Atta (1978); 
and Asfar et al. (1979), has been extended as documented in Chapter 
3. This modified tangential velocity calculation method allows both 
use of non-rectangular panels and movement of the control point 
location away from the panel center. Comparisons of solutions 
obtained using these different formulations have been displayed. It 
is also found that the present modified equation yields the same 
accuracy Cp distribution as compared to Atta's method, but allows 
more flexibility in control point location. None of these methods 
yield satisfactory Cp distributions at the leeward side of an 
ellipsoidal body nose and tail when the flow is at a moderate angle 
of attack. Next, in Asfar et al. (1979), it has been asserted that 
the vortex lattice surface panel method is more accurate than the 
surface source panel method, since the numerical solutions of the 
vortex lattice method matched experimental data for 30° cone- 
cylinder much better than for the source method. In the present 
work, it has been proven that this conclusion was incorrect, being 
due to differences in geometries used in the two methods. Error 
sensitivities of this method to locating the control points off the 
actual body surface, in the plane of the panel, have not been 
investigated in the present work. 

Since the vortex lattice surface panel method has the ability 
to model an arbitrary body, but required a great deal of 
computational effort, a simplified vortex lattice model based on the 
technique used by Mason et al. (1977), using slender body 
transformation of the actual body slope to a bounding rectangular 
box has been investigated as a part of the present research. A 
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brief description of this study has been given by Shu and Kuhlman 
(1981). Due to difficulty in obtaining an adequate expression for 
the local tangential velocity on the rectangular box, this approach 
has been found to not be fruitful. 

Finally, as previously discussed, the main purpose of the 
present study has been to develop theoretical methods which will 
allow the inclusion of fuselage effects into a pair of existing 
subcritlcal wing design programs (Kuhlman and Shu, 1981; Kuhlman, 
1983). Both fuselage programs are now ready to be modified for 
addition into the above mentioned design computer programs. Future 
plans are to accomplish this task. 
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APPENDIXES 



APPENDIX A 

MACSYMA INTEGRAL EXPRESSIONS 


Integral expressions, utilized In the axial singularity method 
of Chapter 2, obtained by using the MACSYMA symbolic manipulation 
language (Bogen et al., 1975) are given below: 


(1) H, = / 


C(x t -5r + zfy 


(x 1 -e) (z 1 2 + x. Z -Zx^ + ? 2 ) 1/2 
z. 4 + (x. Z -2x.? + c 2 )z 2 

1 1 1 Va.d 


(2) Hy = / 


C(x r c ) 4 + z. c Y 


2 2 
z i + x i “V 

z^tz. 2 + x 2 . -2x.$ + ? 2 )' 


(A. 2) 


Ho = / 


t(x r 5) £ + r 


-[z 1 2 (z i 2 +x. 2 -2x.C+C 2 ) 1/2 sinh" 1 (^-^) + (x.+5)z. 2 +x. 3 -x. 2 5] 


z i (z i + x i " 2x i^‘ ) 


(A. 3) 


(4) P x = / 


c(x r 5)‘ + z. 4 r 


( x . -5+2 ) z ] . 2 +2x . 2 -4x . £+2 c 2 


(z i t +x f -2x^+5 ) v C-z^ +('<x i -5x 1 C+3r)z 1 ] 


(A. 4) 


(5) Po = / 


CdJ- 


[(x r c> z + TfF* 


z i 4 +(3x i 2 -3x i 5)z i 2 +2x i 4 -6x. 3 5+6x i 2 5 2 -2x f C 3 

(z. 2 +x i 2 -2x l .5+5 2 ) 1//2 C3z 1 6 +{3x 1 2 -6x i C+35 2 )z 1 4 ] 


(A. 5) 


(6) P, - / 




C(x r 0 2 + z i 2 ] 5/2 


-[2x 1 z i 4 +(4x i 3 -6x i 2 C+3x i 5 2 -5 3 )z i 2 +2x 1 5 -6x i 4 C+6x 1 . 2 5 2 -2x i 2 ? 3 ] 


(z i 2 +x i 2 -2x i ?+c 2 ) 1/2 [3z. 6 +(3x. 2 -6x.c+3c 2 )z/] 




'I 


(A. 6) 


C 3 dC 


(7) P4 = / 2 2-i 5/? 

4 [(x.-c) 2 + z. 2 ] 5/2 


-[2z 1 6 + (6 Xj 2 -6xj 5 + 3e 2 )z t 4 + ( ex^ 4 - 12x^5 + 9x 1 - 2 5 2 -3x 1 5 3 )z 1 2 + 


(z f 2 + Xj Z -2x t 5 + PT^Hz^ + 


2x i - 6 - 6 x^ 5 § + 6x^ 4 C 2 - 2x i - 3 C 3 ] 
(3x^ 2 - 6 x- 5 + 3£ 2 )z i 4 ] 


(A.7 ) 
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APPENDIX B 


SELECTION OF CONTROL POINT LOCATION 
FOR VORTEX LATTICE METHOD 

According to Asfar et al. (1979), the most accurate Cp solutions 
will be obtained when the control points have been located so that 
the self-induced velocity due to the closed loop vortex filament at 
its own control point is a minimum. Searching for the correct 
location of this point, of course, is very tedious and requires a 
lot of calculations. In the present study of axisymmetric bodies, 
all panels have been chosen to be planar quadrilaterals with two 
parallel sides and having the other two sides of equal length, 
except at the nose and tail where panels have been chosen to be of 
equilateral triangular shape. Through a series of numerical 
experiments, as shown below, it has been found that choosing the 
control point location to be the centroid of area of each panel can 
be used to replace the original technique without too much 
difference in the solution accuracy. The following two sections are 
designed to show, first, how to compute the coordinates of the 
centroid of area of a typical panel, and second, the difference in 
location of control points calculated for these two methods. 

B.l Determine the Coordinates of the Centroid of Area 

Consider a typical panel as described above, located in local 
cartesian coordinates as shown in Sketch B.l. 
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Sketch B.1 Nomenclature of a typical panel in 
local cartesian coordinates 


The local coordinates of the centroid of area are given as 


x 


c 


(1 


2a+b 

3(a+b) 


) A 


*c 


= 0 


(B.l ) 
(B.2) 


The global coordinates can be obtained by performing a simple 
transformation. 

B.2 Difference in Coordinates between the Minimum Self-Induced 
Velocity Point and the Centroid of Area 

The following brief table, obtained from numerical experiments, 
shows that the centroid of area may be used in place of the optimal 
control point location, because small differences exist in control 
point location and induced velocity magnitude. 
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Table B.l Comparisons of the minimum self-induced 
velocity point and the centroid of area 


Test No. 

Panel 

Geometry 

Centroid 

of 

Area 


Minimum 

Induced 

Velocity 

Di f ferences 
V\:* W 

a 

b 

l 

x c 

*c 

x m 

*m 

L 

w m 

1 

2. 

4. 

5. 

2.78 

0 

2.87 

0 

1.8 

0.6 

2 

0 . 

4. 

6. 

4. 

0 

4.3 

0 

6.0 

1.2 

3 

0 . 

1.364 

2. 

1.34 

0 

1.45 

0 

5.5 

1.2 

4 

1.364 

1.819 

2. 

1.048 

0 

1.07 

0 

1.1 

0.2 

5 

2.082 

2.229 

2. 

1.011 

0 

1 .04 

0 

1.45 

0.16 

6 

2.229 

2.274 

2. 

1.003 

0 

1.04 

0 

1.83 

0.3 
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APPENDIX C 


CALCULATION OF LOCAL TANGENTIAL VELOCITY 
FOR VORTEX LATTICE METHOD 

When the field piont is on the body surface, one has to account 
for the induced tangential velocity due to the local strength of the 
vortex sheet. The parameters involved in calculating the components 
of the induced tangential velocity in the x and y directions at a 
point P for a rectangular vortex element on the X-Y plane are shown 
in Sketch C.l. With linear interpolation, Kandil et al. (1977) gave 
the expressions for the components of tangential velocity as 
f ol 1 ows : 



Sketch C.l Nomenclature of calculating the induced 
tangential velocity 
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(C.l) 


T vx = *- r l ^ x i +i i~ x ) + r 3 (x-x 1 )]e 1 

V = rr Cr 4 (y i + V y) + (c. 2 ) 

^x>2 


Extension of these two equations to a general, quadrilateral vortex 
element has been given by Atta (1978) as follows: 


aV = (■ 


r l*l + r 2*2 + r 3*3 
AZA 


+ r 4 ^4 


■) xn 


(C.3) 


Note that this expression yields a tangential velocity jump which 
either must be interpreted as being constant on the element, or it 
must be understood to apply only at control points located at the 
average of the panel corner points. 
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Table 1 Comparison of the small a solution by Maruhn with the exact 
solution by Lamb for an ellipsoid of SR=5 along top meridian 

1 ine. 



C p at 

O 

to 

It 

0 


c p 

at a * 10° 


X/L 

MARUHN 

LAMB 

DEVIATION 

MARUHN 

LAMB 

DEVIATION 

.0375 

-.04460 

-.04223 

-.00236 

-.18241 

-.17024 

-.01217 

.0875 

-.12783 

-.12164 

-.00619 

-.22203 

-.19556 

-.02647 

.1375 

-.14472 

-.13736 

-.00736 

-.21531 

-.18466 

-.03065 

.1875 

-.14836 

-.14046 

-.00790 

-.20303 

-.17055 

-.03248 

.2375 

-.14747 

-.13927 

-.00820 

-.19002 

-.15661 

-.03341 

.2875 

-.14453 

-.13615 

-.00839 

-.17712 

-.14322 

-.03390 

.3375 

-.14042 

-.13194 

-.00847 

-.16437 

-.13025 

-.03412 

.3875 

-.13549 

-.12696 

-.00853 

-.15163 

-.11745 

-.03417 

.4375 

-.12984 

-.12130 

-.00854 

-.13866 

-.10457 

-.03409 

.4875 

-.12346 

-.11493 

-.00853 

-.12521 

-.09132 

-.03389 

.5375 

-.11622 

-.10773 

-.00849 

-.11095 

-.07737 

-.03358 

.5875 

-.10789 

-.09948 

-.00841 

-.09546 

-.06232 

-.03314 

.6375 

-.09811 

-.08981 

-.00830 

-.07815 

-.04559 

-.03256 

.6875 

-.08630 

-.07815 

-.00815 

-.05815 

-.02639 

-.03176 

.7375 

-.07145 

-.06353 

-.00792 

-.03408 

-.00340 

-.03068 

.7875 

-.05174 

-.04416 

-.00758 

-.00348 

0.02565 

-.02217 

.8375 

-.02341 

-.01636 

-.00705 

0.03859 

0.06531 

-.02672 

.8875 

0.02292 

0.02900 

-.00608 

0.10393 

0.12643 

-.02250 

.9375 

0.11984 

0.12366 

-.00382 

0.23179 

0.24468 

-.01289 

.9875 

0.53989 

0.53197 

-.00793 

0.71349 

0.67851 

-.03497 
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Table 2 Solution accuracy for axial singularity method for axi symmetric flow 
past an ellipsoid of SR=5. 


Inset 


Config. 

No. 

Axial Flow 

Numerical Arrangement 

Panelling 

a/L*100 

Root 

Mean 

Square 

Error 

Singularity 

Type 

Number of 
Elements 

Panelling 

Scheme 

1 

CS 

10 

Equal 

0 

0.054962 

2 

CS 

10 

Cosine 

1 

0.001758 

3 

CS 

20 

Equal 

0 

0.052254 

4 

CS 

20 

Equal 

0.5 

0.026120 

5 

CS 

20 

Equal 

1 

0.001695 

6 

CS 

20 

Equal 

1.5 

0.025122 

7 

CS 

20 

Equal 

2 

0.040193 

8 

CS 

20 

Cosine 

0 

0.006483 

9 

CS 

20 

Cosine 

't 

1 

0.000089 

10 

CD 

20 

Equal 

1 

0.347660 

11 

LVD 

20 

Equal 

0 

0.052254 

12 

LVD 

20 

Equal 

1 

0.001688 

13 

LVS 

20 

Equal 

1 

0.000503 

14 

LVS 

20 

Cosine 

1 

0.000063 

15 

CS 

20 

Equal 

0.75 

0.012057 

16 

CS 

20 

Equal 

1.75 

0.014300 

17 

CS 

40 

Equal 

1 

0.000062 

18 

CS 

40 

Cosine 

1 

0.000008 

19 

LVS 

10 

Equal 

1 

0.000595 
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DOUBLET STRENGTH SOURCE STRENGTH 





x/L 


(b) Piecewise constant doublet and linearly varying doublet. 


Figure 2. Axial singularity strength distribution along body 
axis for ellipsoidal body, SR ■ 5, a » 0°, equally 
spaced panels and a/L * 0.01. 


94 










Table 3 Comparison of Cp error for various of axial singularity 
types for axi symmetric flow past an ellipsoid of SR=5. 




DEVIATION FROM 

EXACT SOLUTION 

X 10~ 4 

X/L 

CONSTANT 

SOURCE 

LINEARLY 

VARYING 

SOURCE 

CONSTANT 

DOUBLET 

LINEARLY 

VARYING 

DOUBLET 

.0375 

26.209 

-4.542 

-1616.645 

26.209 

.0875 

4.91 

-0.910 

-449.064 

4.909 

.1375 

2.075 

-0.423 

-204.528 

2.075 

.1875 

1.165 

-0.248 

-118.926 

1.165 

.2375 

0.773 

-0.168 

-79.899 

0.774 

.2875 

0.571 

-0.125 

-59.278 

0.570 

.3375 

0.455 

-0.101 

-47.489 

0.455 

.3875 

0.388 

-0.086 

-40.605 

0.388 

.4375 

0.352 

-0.078 

-36.834 

0.352 

.4875 

0.377 

-0.075 

-35.347 

0.337 

.5375 

0.342 

-0.076 

-35.847 

0.342 

.5875 

0.366 

-0.082 

-38.436 

0.366 

.6375 

0.416 

-0.093 

-43.657 

0.415 

.6875 

0.501 

- 0.111 

-52.766 

0.501 

.7375 

0.648 

-0.144 

-68.473 

0.648 

.7875 

0.913 

-0.201 

-97.087 

0.913 

.8375 

1.438 

-0.317 

-155.847 

1.438 

.8875 

2.716 

-0.600 

-306.678 

2.715 

.9375 

7.094 

-1.618 

-921.745 

7.090 

.9875 

70.616 

21.952 

-15422.648 

69.703 
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a/L x 100 


Figure 4 


Effect of inset of axial singularity 
distribution for an ellipsoid of SR 1 
for equally spaced constant source 




Table 4 Solution accuracy for axi symmetric flow past 
ellipsoids, using 20 cosine spaced linearly 
varying axial source panels. 


config. 

inset 

root mean 

No. 

a/Li 

square root 

1 

2 

0.146652 

2 

3 

0.002798 

3 

4 

0.001519 

4 

5 

0.000047 

5 

6 

0.000045 

6 

7 

0.000078 


(a) Solution accuracy for various inset 
values of SR=2 ellipsoid 


config. 

inset 

root mean 

No. 

a/Li 

square error 

1 

0 

0.000330 

2 

0.1 

0.000597 

3 

0.2 

0.000301 

4 

0.25 

0.000004 

5 

0.3 

0.000381 

6 

0.4 

0.001423 

7 

0.6 

0.004844 

8 

0.8 

0.010803 

9 

1.0 

0.019887 


(b) Solution accuracy for various inset 
values for SR=10 ellipsoid 


config. 


inset 

root mean 

No. 

SR 

a/Li 

square error 

1 

2 

6 

0.000045 

2 

3 

2.78 

0.000095 

3 

4 

1.56 

0.000106 

4 

5 

1 

0.000063 

5 

6 

0.7 

0.000005 

6 

7 

0.51 

0.000022 

7 

8 

0.4 

0.000061 

8 

9 

0.3 

0.000068 

9 

10 

0.25 

0.000004 


(c) Solution accuracy for various SR 
ellipsoids, using the suggested 
Inset value calculated by equation 
(2.69) 



EXACT SOLUTION 






Figure 6. Pressure coefficient distribution along meridian 
line for ellipsoid of SR - 5, at a - 0°, cosine 
spaced linearly varying axial source with 
a/L ■ 0.01 (optimal panelling). 
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ERROR 


5 



Figure 7. Error distribution along meridian line for an ellipsoid 
of SR m 5, at c l ■ 0°, cosine spaced linearly varying 
axial source with a/L - 0.01 (optimal panelling). 
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Table 5 Solution accuracy for axial singularity method for inclined flow 
with 5° angle of attack past an ellipsoid of SR=5; error 
calculated along 0 s 33.75° meridian line on body surface. 


CROSS FLOW-NUMERICAL ARRANGEMENT 


CONFIG. 

NO. 

SINGULARITY 

TYPE 

ELEMENT 

NO. 

PANELLING 

SCHEME 

( a/L )*1 00 

ROOT MEAN 
SQUARE ERROR 

1 

CD* 

20 

EQUAL 

0 

0.019056 

2 

CD* 

20 

EQUAL 

1 

0.010986 

3 

CD* 

20 

COSINE 

0 

0.002012 

4 

CD* 

20 

COSINE 

1 

0.000117 

5 

LVD* 

20 

EQUAL 

0 

0.005749 

6 

LVD* 

20 

EQUAL 

1 

0.000129 

7 

LVD* 

20 

COSINE 

0 

0.011609 

8 

LVD* 

20 

COSINE 

1 

0.000088 

9 

LVD + 

20 

COSINE 

1 

0.000062 

* - CONFIGURATIONS BASED 

ON FIXED 

CONDITIONS IN 

AXIAL FLOW: 

20 DISCRETE 


ELEMENTS, COSINE SPACING, CONSTANT SOURCE WITH a/L = 0.01. 

+ - AXIAL FLOW CHANGED TO 20 DISCRETE ELEMENTS, COSINE SPACING, LINEARLY 
VARYING SOURCE WITH a/L = 0.01, THE OPTIMAL NUMERICAL ARRANGEMENT FOR 
INCLINED FLOW. 
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DOUBLET STRENGTH SOURCE STRENGTH 



(a) Axial flow: 20 discrete elements, cosine spaced linearly 

varying source with a/L ■ 0.01. 



x/L 

(b) Cross flow: 20 discrete elements, cosine spaced linearly 

varying doublet with a/L » 0.01. 


Figure 8. Axial singularity strength distribution along body 
axis, for an ellipsoid of SR * 5, at a » 30°, for 
optimal panelling. 
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EXACT SOLUTION 


Figure 9 


AXIAL SINGULARITY METHOD 


Pressure coefficient distribution from axial 
singularity method for an ellipsoid of SR = 5, 
at a » 30°, for optimal panelling, various meridian 
lines on body;- 9-0°, 9 => 30°, 9 - 60° and 9 - 90 ° 




EXACT SOLUTION 


AXIAL SINGULARITY METHOD 


Figure 12. Pressure coefficient distribution from axial 

singularity method along top meridian line (9 » 0°) 
on body surface, for an ellipsoid of SR » 5, at 
o*0,oi a 10,oi a 20 and ct * 30° for the optimal 
panelling. 





Table 6 Comparison of Cp from vortex lattice method for various of 

tangential velocity calculation methods for axisymmetric flow 
past an ellipsoid of SR=5, using 8*30 panels, cosine spacing 
axially, and control point located at the average of four 
corner coordinates of each panel. 



Cp 


Cp NUMERICAL 


x/L 

EXACT 

KANDIL 

(1977) 

ATTA 

(1978) 

MODIFIDIED 
(equ. 3.20) 

.00137 

.864390 

.9329 

.8387 

.8387 

.01770 

.269314 

.2868 

.2607 

.2607 

.05511 

.026309 

.0276 

.0180 

.0180 

.11196 

-.057666 

-.0629 

-.0670 

-.0670 

.18577 

-.093193 

-.1016 

-.1035 

-.1035 

.27332 

-.110248 

-.1203 

-.1211 

-.1211 

.37077 

-.118535 

-.1295 

-.1297 

-.1297 

.47387 

-.121615 

-.1329 

-.1329 

-.1329 

.57811 

-.120616 

-.1318 

-.1319 

-.1319 

.67894 

-.115183 

-.1258 

-.1262 

-.1262 

.77195 

-.103200 

-.1126 

-.1138 

-.1138 

.85307 

-.078821 

-.0859 

-.0887 

-.0887 

.91876 

-.025376 

-.0280 

-.0342 

-.0342 

.96616 

.113794 

.1211 

.1057 

.1057 

.99317 

.539069 

.5755 

.5286 

.5286 
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Table 7 Comparison of Cp from vortex lattice method for various of 

tangential velocity calculation methods for axlsymmetrlc flow 
past an ellipsoid of SR=5, using 8*30 panels, cosine spacing 
axially, and control point located at the centroid of area of 
each panel . 





Cp NUMERICAL 





SCHEME 1* 


SCHEME 

2 + 

X/L 

Cp 

EXACT 

KANDIL 

(1977) 

ATTA 

(1978) 

MODIFIED 
(equ 3.20) 

ATTA 

(1978) 

MODIFIED 
(equ. 3.20) 

.00183 

.825989 

.9018 

.8446 

.8446 



.01814 

.262712 

.2844 

.2680 

.2601 

.2623 

.2704 

.05551 

.025187 

.0320 

.0260 

.0231 

.0234 

.0264 

.11231 

-.057943 

-.0560 

-.0585 

-.0598 

-.0597 

-.0584 

.18606 

-.093277 

-.0936 

-.0947 

-.0953 

-.0953 

-.0947 

.27352 

-.110274 

-.1118 

-.1122 

-.1125 

-.1124 

-.1122 

.37089 

-.118541 

-.1207 

-.1208 

-.1209 

-.1208 

-.1208 

.47389 

-.121615 

-.1240 

-.1240 

-.1240 

-.1240 

-.1240 

.57804 

-.120619 

-.1229 

-.1229 

-.1230 

-.1230 

-.1229 

.67877 

-.115199 

-.1171 

-.1173 

-.1175 

-.1174 

-.1173 

.77170 

-.103247 

-.1043 

-.1050 

-.1054 

-.1053 

-.1050 

.85275 

-.078971 

-.0783 

-.0800 

-.0809 

-.0808 

-.0800 

.91838 

-.025914 

-.0220 

-.0259 

-.0278 

-.0276 

-.0257 

.96572 

.111215 

.1229 

.1133 

.1086 

.1094 

.1141 

.99272 

.521222 

.5696 

.5385 

.5248 

.5397 

.5533 

* - CONTROL POINT 

LOCATED AT THE 

CENTROID 

OF AREA OF 

EACH PANEL. 


+ - CONTROL POINT LOCATED AT THE CENTROID OF AREA OF EACH PANEL, 
EXCEPT TRIANGULAR PANEL AT THE NOSE AND TAIL, WHERE THE CONTROL 
POINT LOCATED AT 3/4 CHORD LENGTH AWAY FROM THE APEX. 
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Table 8 Comparison of Cp from vortex lattice method for various of 
tangential velocity calculation methods for inclined flow, 
a=20°, past an ellipsoid of $R s 5, using 12*20 panels, 
cosine spacing axially, and control point located at the 
average of four corner coordinates of each panel. 



e 

= 7.5° 


0 3 

82.5° 


x/L 

EXACT 

ATTA 

(1978) 

MODIFIED 
(equ. 3.20) 

EXACT 

ATTA 

(1978) 

MODIFIED 
(equ. 3.20) 

.003080 

-.100442 

-.2543 

-.2451 

.275594 

.3048 

.3077 

.015320 

-.390676 

-.4691 

-.4609 

-.107441 

-.1122 

-.1101 

.039485 

-.386188 

-.4482 

-.4426 

-.288417 

-.2935 

-.2922 

.074995 

-.316861 

-.3756 

-.3717 

-.359549 

-.3643 

-.3636 

.12097 

-.250026 

-.3066 

-.3038 

-.389816 

-.3946 

-.3941 

.176275 

-.192597 

-.2472 

-.2452 

-.403316 

-.4081 

-.4078 

.239555 

-.142936 

-.1956 

-.1942 

-.408920 

-.4137 

-.4135 

.309250 

-.098582 

-.1493 

-.1483 

-.410256 

-.4150 

-.4148 

.383635 

-.057401 

-.1062 

-.1056 

-.408901 

-.4135 

-.4134 

.460890 

-.017586 

-.0643 

-.0641 

-.405509 

-.4099 

-.4089 

.539110 

.022528 

-.0221 

-.0222 

-.400228 

-.4043 

-.4043 

.616365 

.064701 

.0224 

.0291 

-.392826 

-.3965 

-.3966 

.690750 

.111064 

.0714 

.0705 

-.382655 

-.3859 

-.3860 

.760445 

.164558 

.1279 

.1267 

-.368438 

-.3711 

-.3713 

.823725 

.229688 

.1967 

.1950 

-.347721 

-.3497 

-.3499 

.879030 

.313947 

.2854 

.2833 

-.315568 

-.3167 

-.3168 

.925005 

.430552 

.4078 

.4052 

-.261150 

-.2612 

-.2613 

.960515 

.602127 

.5874 

.5843 

.158303 

-.1571 

-.1570 

.194685 

.846020 

.8420 

.8393 

.056286 

.0604 

.0611 

.996920 

.986992 

.9907 

.9913 

.418758 

.4687 

.4705 
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Table 9 Comparison of Cp from vortex lattice method along 0 = 7.5° 

meridian line for various of tangential velocity calculation 
methods for inclined flow, a = 20°, past an ellipsoid of 
SR=5, using 12*20 panels, cosine spacing axially, and control 
point located at the centroid of area of each panel. 


Cp 

x/L EXACT 


004105 

.171617 

016315 

-.401517 

040435 

-.384262 

07587 

-.315345 

12175 

-.249073 

176945 

-.192009 

240095 

-.142561 

30964 

-.098351 

383875 

-.057274 

46097 

-.017545 

53903 

-.022487 

616125 

.064563 

69036 

.110795 

759905 

.164089 

823055 

.228870 

87825 

.312452 

92413 

.427582 

959565 

.595590 

983685 

.831918 

995895 

.992848 


SCHEME 1* 


ATTA 

(1978) 

MODIFIED 
(equ. 3.20) 

.0539 

.0363 

-.3916 

-.4074 

-.4131 

-.4166 

-.3506 

-.3510 

-.2858 

-.2852 

-.2287 

-.2277 

-.1785 

-.1776 

-.1332 

-.1324 

-.0908 

-.0903 

-.0496 

-.0495 

-.0080 

-.0082 

.0359 

.0352 

.0841 

.0829 

.1389 

.1379 

.2075 

.2046 

.2947 

.2905 

.4146 

.4086 

.5899 

.5810 

.8364 

.8252 

.9903 

.9932 


SCHEME 2 + 


ATTA 

(1978) 

MODIFIED - 
(equ. 3.20) 

-.3412 

-.3639 

-.4092 

-.4131 

-.3496 

-.3500 

-.2854 

-.2848 

-.2285 

-.2275 

-.1783 

-.1774 

-.1331 

-.1323 

-.0907 

-.0902 

-.0496 

-.0494 

-.0079 

-.0081 

.0359 

.0353 

.0842 

.0830 

.1399 

.1380 

.2076 

.2047 

.2949 

.2907 

.4150 

.4089 

.5905 

.5817 

.8361 

.8250 


* - CONTROL POINT LOCATED AT THE CENTROID OF AREA OF EACH PANEL. 

+ - CONTROL POINT LOCATED AT THE CENTROID OF AREA OF EACH PANEL 
EXCEPT TRIANGULAR PANEL AT THE NOSE AND TAIL, WHERE THE 
CONTROL POINT LOCATED AT 3/4 CHORD LENGTH AWAY FROM THE APEX. 
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Table 10 Comparison of Cp from vortex lattice method along 0 ■ 82.5° 
meridian line for various of tangential velocity calculation 
methods for inclined flow, a = 20°, past an ellipsoid of 
SR=5, using 12*20 panels, cosine spacing axially, and control 
point located at the centroid of area of each panel. 




SCHEME 

1* 

SCHEME 

2 + 

x/L 

Cp 

EXACT 

ATTA 

(1978) 

MODIFIED 
(equ. 3.20) 

ATTA 

(1978) 

"""MODI? 1 ED 
(equ. 3.20) 

004105 

.214953 

.1146 

.2799 



016315 

-.122094 

-.1164 

-.1162 

-.0979 

-.1002 

040435 

-.291783 

-.2878 

-.2896 

-.2859 

-.2878 

07587 

-.360491 

-.3545 

-.3555 

-.3539 

-.3550 

12175 

-.390117 

-.3827 

-.3832 

-.3824 

-.3830 

176945 

-.403413 

-.3951 

-.3954 

-.3950 

-.3952 

240095 

-.408945 

-.4001 

-.4002 

-.4000 

-.4000 

30964 

-.410255 

-.4010 

-.4010 

-.4009 

-.4009 

383875 

-.408893 

-.3993 

-.3992 

-.3992 

-.3992 

46097 

-.405505 

-.3956 

-.3956 

-.3956 

-.3956 

53903 

-.400234 

-.3901 

-.3902 

-.3901 

-.3901 

616125 

-.392853 

-.3826 

-.3827 

-.3825 

-.3826 

69036 

-.382720 

-.3724 

-.3726 

-.3723 

-.3725 

759905 

-.368574 

-.3582 

-.3586 

-.3581 

-.3585 

823055 

-.348004 

-.3377 

-.3385 

-.3375 

-.3383 

87825 

-.316190 

-.3063 

-.3075 

-.3060 

-.3072 

92413 

-.262683 

-.2537 

-.2555 

-.2532 

-.2550 

959565 

-.162783 

-.1558 

-.1582 

-.1543 

-.1569 

983685 

.040291 

.0453 

.0460 

.0571 

.0563 

995895 

.368258 

.2520 

.4155 




* - CONTROL POINT LOCATED AT THE CENTROID OF AREA OF EACH PANEL. 

+ - CONTROL POINT LOCATED AT THE CENTROID OF AREA OF EACH PANEL 
EXCEPT TRIANGULAR PANEL AT THE NOSE AND TAIL, WHERE THE 
CONTROL POINT LOCATED AT 3/4 CHORD LENGTH AWAY FROM THE APEX. 
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.0 


Figure 14. Pressure coefficient distributions from vortex lattice method 
along various meridian lines for ellipsoidal body, SR a 5, 
a * 20° , using modified tangential velocity calculation 
method . 



115 


EXACT SOLUTION 


VORTEX LATTICE SURFACE PANEL METHOD; 12 x 20 PANELS 



6 (DECREES) 

(d) x/L - 0.53903. 


0 (DECREES) 

(e) x/L - 0.69036. 


0 (DECREES) 

(f) x/L » 0.92413. 


Figure 15. Pressure coefficient distributions from vortex lattice method along circumferential 
lines at various x/L values for ellipsoidal body, SR = 5, a « 20°, using modified 
tangential velocity calculation method. 



STRENGTH 


CONSTANT SOURCE DISTRIBUTION 


LINEARLY VARYING SOURCE DISTRIBUTION 

A 

/ I 
/ I 


1 i /n- 

0.4 I-.,/ 0.6 

1 / 

x/L * 


(a) Source strength distribution along body axis. 



x 3/2 + (5y) 2 - 1 


(b) Pressure distribution along meridian line. 


Axial singularity strength, and C^ distribution, for 

a SR » 5, modified ellipsoid given by equation (4.1), 
using 20 cosine spaced, constant and linearly varying 
axial sources. 



Figure 17. Pressure coefficient distribution for a SR ■ 5, modified 

ellipsoid given by equation (4.1), using axial singularity 
method and vortex lattice surface panel method. 
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STRENGTH 


0.0 0.2 0.4 0.6 0.8 1.0 

x/L 

0.0 0.2 0.4 0.6 0.8 1.0 

x/L 

(a) Strength distribution 
using constant source 
panels . 

(b) Strength distribution 
using linearly varying 
source panels. 



•0.8 20 EQUALLY SPACED LINEARLY VARYING SOURCE PANELS 

■1.0* 20 EQUALLY SPACED CONSTANT SOURCE PANELS 

(c) Pressure coefficient distribution along meridian line. 

Figure 18. Axial singularity strength, and C distribution, for 

P 

axisymmetric flow past a SR » 4, Ogival body of 
revolution, using 20 equal spaced, constant and 
linearly varying axial sources. 
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AXIAL SINGULARITY METHOD; 20 EQUAL SPACED LINEARLY 

VARYING SOURCE PANELS 

VORTEX LATTICE SURFACE METHOD; 8 * 30; COSINE 

AXIAL SPACING 

VORTEX LATTICE SURFACE METHOD; 8 * 30; EQUAL 

AXIAL SPACING 

Figure 19. Pressure coefficient distribution for axisymmetric 

flow past a SR ■ 4, Ogival body of revolution, using 
axial singularity method and vortex lattice surface 
panel method. 


STRENGTH 




x/L x/L 


(a) Strength distribution (b) Strength distribution 
using constant source using linearly varying 

panels. source panels. 



Figure 20. Pressure coefficient distribution along meridian line 
for axisymmetric flow past a cylindrical body with 
modified ellipsoidal nose and tail, using axial 
singularity method, and vortex lattice surface 
panel method. 
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EXPERIMENTAL DATA BY FOX (1971); M - 0.403 

— AXIAL SINGULARITY METHOD; 10-1-10 LINEARLY 
VARYING SOURCE PANELS 

VORTEX LATTICE SURFACE PANEL METHOD; 

8 x 30 COSINE AXIAL SPACING 


COMPUTATIONAL MODEL FOR AXIAL 
SINGULARITY METHOD AND VORTEX 
LATTICE SURFACE PANEL METHOD - 


Figure 21. Pressure coefficient distribution along meridian line for 
axisymmetric flow past a cylindrical body with Ogival nose 
using axial singularity method, and vortex lattice surface 
panel method. 
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Figure 22. Pressure coefficient distribution along meridian line for axisymmetric 
flow past a cylindrical body with SR ** 2 ellipsoidal nose, using axial 
singularity method, and vortex lattice surface panel method. 




AXIAL SINGULARITY METHOD; AXIAL FLOW: 10-3-10 PANELS 

LINEARLY VARYING SOURCE; CROSSFLOW: 10- 3-10 PANELS 

LINEARLY VARYING DOUBLET 

VORTEX LATTICE SURFACE PANEL METHOD; 10 * 32 COSINE 

AXIAL SPACING; 8-9° 


Figure 23. Pressure coefficient distribution along top meridian line for 
inclined flow a ■ 6.08° past a cylindrical body with SR ■ 2 
ellipsoidal nose, using axial singularity method, and vortex 
lattice surface panel method. 
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EXPERIMENTAL DATA BY JOHNSON (1963); 
L/R =■ 4.8; M - 0.18 
SURFACE SOURCE METHOD OF HESS (1967) 
SEMI-INFINITE BODY 

SURFACE VORTEX LATTICE METHOD; 

L/R - 4.8 (PRESENT) 

SURFACE VORTEX LATTICE METHOD: 

L/R - 9.2 (PRESENT) 
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Figure 24. Pressure coefficient distribution along meridian line 
for axisymmetric flow past a 30° cone-cylinder for 
vortex lattice surface panel method. 
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EXPERIMENTAL DATA BY JOHNSON (1963): 

M » 0,18 

SURFACE SOURCE METHOD OF HESS (1967); 

SEMI-INFINITE BODY 

SURFACE VORTEX LATTICE METHOD; 

L/R » 4.8 (PRESENT) 
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Pressure coefficient distribution along 9-90 
meridian line for flow at a * 20° past a 30° 
cone-cylinder for vortex lattice surface panel 
method. 
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